iteration to convergence

I am very very new to Matlab and this will become quite apparent, but I am working on a problem where the variable I'm solving for is also in the equation. I was told this could be accomplished in Matlab.
v=(((2*Ke)/deltaZ)*log(Ar))+ sqrt((alpha+v^2)/2)
PS. this is a Hydrogeology equation to determine velocity using heat as a tracer
%Constants
lambdaE=2.8; %Effecitve Thermal Conductivity
p=2.7; %Density of Fluid and Sediment
c=4.2; %Heat Capacity of Fluid and Sediment
Ke=[lambdaE/(p*c)]; %Effective Thermal Diffusivity
f=2; %Frequency
P=(1/f); %Period
v=.3; %Fluid Velocity (initial)
deltaPhi = 0.004; %Phase shift between shallow and deep points [Measured in the Lab]
Ar = 0.1; %Amp ratio of shallow and deep points[Measured in the Lab]
deltaZ=.2; %Distance between shallow and deep points
alpha=sqrt(v^4+(8*pi*Ke/P)^2);
Amp Method
%Ar=exp((deltaZ/(2*Ke))*(v-sqrt((alpha + v^2)/2)))
%VAr
v=(((2*Ke)/deltaZ)*log(Ar))+ sqrt((alpha+v^2)/2)
I really hope someone out there can take pity and give me a brief walk through how this can be solved.
-Ethan

 Accepted Answer

bym
bym on 11 Sep 2011
I think it converges pretty quick; see below. My tweaks are commented with % *comment *
clc;clear
%%Constants
lambdaE=2.8; %Effecitve Thermal Conductivity
p=2.7; %Density of Fluid and Sediment
c=4.2; %Heat Capacity of Fluid and Sediment
Ke=lambdaE/(p*c); %Effective Thermal Diffusivity
f=2; %Frequency
P=(1/f); %Period
deltaPhi = 0.004; %Phase shift between shallow and deep points [Measured in the Lab]
Ar = 0.1; %Amp ratio of shallow and deep points[Measured in the Lab]
deltaZ=.2; %Distance between shallow and deep points
Amp Method
%Ar=exp((deltaZ/(2*Ke))*(v-sqrt((alpha + v^2)/2)))
%VAr
tol = 1;
v=zeros(50,1); % *changed v to vector*
alpha = v; % *changed alpha to vector*
v(1)=.3; %Fluid Velocity (initial)
alpha(1)=sqrt(v(1)^4+(8*pi*Ke/P)^2);
c = 0; % *counter*
while tol>.00001 % *desired accuracy*
c=c+1;
v(c+1)=(((2*Ke)/deltaZ)*log(Ar))+ ...
sqrt((alpha(c)+v(c).^2)/2); % * note changes to alpha & c
alpha(c+1)=sqrt(v(c+1)^4+(8*pi*Ke/P)^2); % *added*
tol = abs(v(c)-v(c+1));
end

1 Comment

But what if soluion diverges then the while loop will run till infinity how to break that

Sign in to comment.

More Answers (1)

The solutions are
(4*Ke*ln(Ar) +/- (8*Ke^2*ln(Ar)^2+alpha*deltaZ^2)^(1/2)) / deltaZ
That is, can be done as a simple quadratic.

4 Comments

x y
x y on 8 Sep 2011
Wow! Well I take my hat off to you for your mad math skills, but its more of a "how to do iteration" question. Thanks though, that neat to see the equation reworked like that!
My "mad math skills" just consisted of pasting it in to a symbolic mathematics program ;-)
while true
new_v=(((2*Ke)/deltaZ)*log(Ar))+ sqrt((alpha+v^2)/2);
if abs(new_v - v) < 0.0001 %or as appropriate for tolerance
break; %we converged
end
v = new_v;
end
x y
x y on 11 Sep 2011
Sorry, didn't notice you'd updated the comment. I understand what you did and thats a brilliant idea, but it doesn't work? doesn't seem to iterate nor converge.

Sign in to comment.

Categories

Asked:

x y
on 8 Sep 2011

Commented:

on 5 Sep 2018

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!