# for loop with stored variables

5 views (last 30 days)

Show older comments

Hi i have written a for loop for several variables that i have defined in other coding. I cannot seem to get matlab to recognise these vectors that i have already defined, I am relatively new to matlab and i am unsure if i have set up the for loop correctly. Any help would be great

this defines daisyworld 1

function dydt =daisyworld1(t,y,Lt,L)

Lint=interp1(Lt,L,t); %interpolate the data set (Lt,L) at time t

dydt = [0;0];

A=((1-y(1)-y(2))*0.5)+(y(1)*0.25)+(y(2)*0.75); %albedo

S=917; %constant solar energy

Z=5.67*10^(-8); %Stefan-Boltzmann constant;

Te=((((S*Lint)/Z)*(1-A)).^(1/4))-273; %plantary temperature

B= 1-0.003265*((22.5-Te).^2); %beta value,

g=0.3; %death term

dydt(1)=y(1)*(( (1-y(1)-y(2)) *B)-g); %black daisy formula

dydt(2)=y(2)*(( (1-y(1)-y(2)) *B)-g); %white daisy formula

I then solve it using ode45

clear; % Remove stored variables daisyworldode45

Lt=linspace(0,1000,25); %generate t for L

L=Lt/500+ 0.4; %luminosity function want to keep within a range of 0 and 2 so change depending on tspan

tspan=[0 1000]; %solve for values of t

IC = [0.2 0.3]; %initial conditions of daisy percentage, black then white

[T, Y]=ode45(@(t,y)daisyworld1(t,y,Lt,L),tspan,IC); % solves equation, need capital T and Y

I am trying to extract these values of T and Y, and use the vector T to find values for L and then plot these against Te

N=7269;

Te= zeros(1,N);

for T=0:N

for Y=0:N

Te(T,Y)=((917*(T/500+0.4)/(5.67*10^(-8))*(1-(((1-y(1)-y(2))*0.5)+(y(1)*0.25)+(y(2)*0.75))).^(1/4)))-273;

end

end

plot (T,Y);

7269 is the number of y values that ode45 uses i keep getting the error code

Undefined function 'y' for input arguments of type 'double'.

Error in daisyworldode45 (line 23)

Te(T,Y)=((917*(T/500+0.4)/(5.67*10^(-8))*(1-(((1-y(1)-y(2))*0.5)+(y(1)*0.25)+(y(2)*0.75))).^(1/4)))-273;

##### 3 Comments

Jan
on 14 Mar 2013

### Accepted Answer

ChristianW
on 14 Mar 2013

There are several ways. I would make for the Plantary Temperature (Te) a new function, and vectorize Albedo (A). Like this:

function test

clear % Remove stored variables daisyworldode45

tspan = [0 1000]; %solve for values of t

L = @(t) t/range(tspan)*2 + 0.4; %luminosity function want to keep within a range of 0 and 2 so change depending on tspan

IC = [0.2 0.3]; %initial conditions of daisy percentage, black then white

[T,Y] = ode45(@(t,y)daisyworld1(t,y,L),tspan,IC); % solves equation, need capital T and Y

LT = L(T);

Te = Temp(LT,Y);

subplot(211), plot(T,Y), xlabel('T, Time'), legend('Y_1','Y_2')

subplot(212), plot(LT,Te)

xlabel('L(T), Luminosity'), ylabel('Te, Plantary Temperature')

function dydt = daisyworld1(t,y,L)

dydt = [0;0];

Lt = L(t); % calculate L at time t

Te = Temp(Lt,y);

B = 1-0.003265*((22.5-Te).^2); % beta value

g = 0.3; % death term

dydt(1) = y(1)*(( (1-y(1)-y(2)) *B)-g); % black daisy formula

dydt(2) = y(2)*(( (1-y(1)-y(2)) *B)-g); % white daisy formula

function Te = Temp(Lt,y)

if isvector(y), y = y(:)'; end % make row vector

A = ((1-y(:,1)-y(:,2))*0.5)+(y(:,1)*0.25)+(y(:,2)*0.75); % albedo

S = 917; % constant solar energy

Z = 5.67*10^(-8); % Stefan-Boltzmann constant;

Te = ((((S*Lt)/Z).*(1-A)).^(1/4))-273; % plantary temperature

### More Answers (0)

### See Also

### Community Treasure Hunt

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

Start Hunting!