matlab program help looping over t and the indices of a matrix to fill in entries
1 view (last 30 days)
Show older comments
so this program seems just to iterate over the first row i need it to iterate over the matrix, and since i couldnt get a true summation to work im adding matrices ad each iteration of t until 1500 so i can get a reasonable approximation here is the code
nx=21; %Number of steps in space(x)
ny=26; %Number of steps in space(y)
niter=1500; %Number of iterations
dx=.1; %Width of space step(x)
dy=.1; %Width of space step(y)
x=0:dx:2; %Range of x(0,2) and specifying the grid points
y=0:dy:2.5; %Range of y(0,2.5) and specifying the grid points
t=1; %Initialize t
%Initial Conditions
p=zeros(ny,nx); %Preallocating p
pn=zeros(ny,nx); %Preallocating pn
ps=zeros(ny,nx); %Preallocating ps
%Analytical solution matrix formation
%Explicit iterative scheme with indices and t
for j=1:nx
for k=1:ny
for t=1:niter
p(k,j)=(((4*1.35)/pi)*(sin(((2*t)-1)*pi*(j-1)/20) ...
*sinh(((2*t)-1)*pi*(k-1)/20)/sinh(((2*t)-1)*pi*(20/25))));
ps(k,j)=pn(k,j)+p(k,j);
pn(k,j)=ps(k,j);
end
end
end
%Plotting the solution
contour(pn)
when i look at the variable pn or ps it has NaN values for everything except the first row
Answers (4)
Roger Stafford
on 26 Nov 2014
When you get to the third row where k = 3 and where t has ascended to 1500, the value of
sinh(((2*t)-1)*pi*(k-1)/20)
in the numerator and
sinh((2*t-1)*pi*20/25)
in the denominator have both overflowed to plus infinity and their quotient would therefore be a NaN. To avoid this, you should use the identity
sinh(A)/sinh(B) = (exp(A-B)-exp(-A-B))/(1-exp(-2*B))
so that applying this to
sinh((2*t-1)*pi*(k-1)/20)/sinh((2*t-1)*pi*20/25)
will no longer produce a NaN. That will get you down to about the 18th row. After that even this quantity produces an overflow and you will thereafter get infinities in p, ps, and pn. No identity can prevent infinities in these last eight rows, but at least there will be no NaNs.
1 Comment
Roger Stafford
on 26 Nov 2014
When I said "No identity can prevent infinities in these last eight rows, but at least there will be no NaNs" I forgotten that the sine factor will be changing the signs of the infinity quantities and this will also produce a NaN from the operation infinity minus infinity.
Image Analyst
on 26 Nov 2014
Standard debugging procedure is to break up that ungodly long equation into smaller individual terms. Then check each term to see if it's valid. Step through the code with the debugger examining terms to see when/where/why one or more of the terms is nan.
0 Comments
christopher zangler-scaduto
on 26 Nov 2014
1 Comment
Roger Stafford
on 27 Nov 2014
As I have already indicated to you, Christopher, by the time you have reached the 18th row or thereabouts, the ratio of the two hyperbolic sines will have exceeded the maximum value allowed for double precision floating point numbers. It is a fundamental difficulty of your algorithm. There is no remedy for this other than to use higher precision numbers, which can be done by using symbolic numbers of the Symbolic Toolbox, or by using for example John D'Errico's high precision package in the file exchange.
See Also
Categories
Find more on Matrix Indexing in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!