Help me, Error message: "?? Subscript indices must either be real positive integers or logicals"

Dear Experts, I'm new in using matlab and I'm not using any programming software before. I wanted to use gaver-stehfest algorithm in inverting a laplace equation.

I made these function codes in three m-files:

Function 1:

function [f,y]=prnfr(s) ome=1; lam=1; f=(ome(1-ome)*s+lam)/((1-ome)*s+lam); y=besselk(0,(sqrt(s*f)))/(s*sqrt(s*f)*besselk(1,(sqrt(s*f))));

Function 2: % % ilt=gavsteh(funname,t,L) % % funname The name of the function to be transformed. % t The transform argument (usually a snapshot of time). % ilt The value of the inverse transform % L number of coefficient ---> depends on computer word length used % (examples: L=8, 10, 12, 14, 16, so on..) % % Wahyu Srigutomo % Physics Department, Bandung Institute of Tech., Indonesia, 2006 % Numerical Inverse Laplace Transform using Gaver-Stehfest method % %Refferences: % 1. Villinger, H., 1985, Solving cylindrical geothermal problems using % Gaver-Stehfest inverse Laplace transform, Geophysics, vol. 50 no. 10 p. % 1581-1587 % 2. Stehfest, H., 1970, Algorithm 368: Numerical inversion of Laplace transform, % Communication of the ACM, vol. 13 no. 1 p. 47-49 % % Simple (and yet rush) examples included in functions fun1 and fun2 with % their comparisons to the exact value (use testgs.m to run the examples) function ilt=gavsteh(funname,t,L)

nn2 = L/2; nn21= nn2+1;

for n = 1:L z = 0.0; for k = floor( ( n + 1 ) / 2 ):min(n,nn2) z = z + ((k^nn2)*factorial(2*k))/ ... (factorial(nn2-k)*factorial(k)*factorial(k-1)* ... factorial(n-k)*factorial(2*k - n)); end v(n)=(-1)^(n+nn2)*z; end

sum = 0.0; ln2_on_t = log(2.0) / t; for n = 1:L p = n * ln2_on_t; sum = sum + v(n) * feval(funname,p); end ilt = sum * ln2_on_t;

Function 3: % Applying Gaver-stehvest to the laplace equation

function y=gs(L); L=input('L = ');

for l=1:L t(l)=l*1; F(l)=gavsteh('prnfr',t(l),L);

end Result1=[t' F'] plot(t,F)

I run gs.m file and I always get this error message: >> gs L = 20 ??? Subscript indices must either be real positive integers or logicals.

Error in ==> prnfr at 4 f=(ome(1-ome)*s+lam)/((1-ome)*s+lam);

Error in ==> gavsteh at 42 sum = sum + v(n) * feval(funname,p);

Error in ==> gs at 9 F(l)=gavsteh('prnfr',t(l),L);

Can anybody help please?

[Walter Roberson]

Please edit the question, select the code, and click on the 'Code {}' button. That will format the code so that it is readable by others.

Answer 1

[Walter Roberson]

You have ome=1; f=(ome(1-ome)*s+lam)/((1-ome)*s+lam); Did you mean f=(ome*(1-ome)*s+lam)/((1-ome)*s+lam); ?

Answer 2

[Albaihaqi Albaihaqi]

My bad, thanks..., yes i mean, ome*(1-ome) Why Matlab doesntmake such a warning for this..

Sorry I'm new using Matlab

Thank you for your answer and advice on using 'code{}' button... next time i will