How can I create a script that solves a function for temperatures over a length?

Below is the function code, followed by the script and obtained error. I'm looking to create a script that generates a vector for 'Ts' and 'Tg' that contains different temperatures along a given z length (say 0 - 0.42m). This ideally would be kept in a function code as this will be built on to give concentration profiles using rate of rxn equations which depend on temperature values.

FUNCTION CODE:

function temp = temperature(z, T)
Ts = T(1);
Tg = T(2);
% temperature profiles
Ts1 = ((-3508440.83).*(z.^5))+((3221859.34).*(z.^4))-((1089760.82).*(z.^3))+((165361.13).*(z.^2))-((12015.96).*z)+(1179.21);
Ts2 = ((-1906981.2).*(z.^3))+((2073664.15).*(z.^2))-((754331.65).*(z))+(92096.47);
Tg1 = ((2796881.8).*(z.^5))-((2553240.18).*(z.^4))+((852324.51).*(z.^3))-((124780.34).*(z.^2))+((6092.92).*(z))+(889.33);
Tg2 = ((-1567143.11).*(z.^3))+((1751141.7).*(z.^2))-((653704.46).*(z))+(81812.04);
% Solid Temperature Loop
if z >=0 && z<=0.32
    Ts = Ts1;
elseif z<=0.407
    Ts = Ts2;
end
% Gas Temperature Loop
if z >=0 && z<=0.32
    Tg = Tg1;
elseif z<=0.407
    Tg = Tg2;
end
temp = [Ts, Tg];
end

SCRIPT

z0 = [0 0.415];
T0 = [0 0];
temp = feval(@temperature,z0,T0);
Ts = temp(1)
Tg = temp(2)

ERROR obtained:

Operands to the logical and (&&) and or (||) operators must be convertible to logical scalar values.

Error in temperature (line 13)

if z >=0 && z<=0.32

Error in Tempprofile (line 4)

temp = feval(@temperature,z0,T0);

Answer 1

[Torsten]

z is a vector.

If you use short-circuit AND (&&) , each expression in the if condition must evaluate to a scalar. This is not the case here. So using && will result in an error.

By the way:

What do you want to test with

if z>= 0 && z <=0.32

for z being a vector ?

[Jenna Broders]

Thanks Torsten, so how can I evaluate the loop using a vector instead?I I've done it before with ODEs in a similar manner and it worked ok. Or would oI have to create a vector column first with 1 column and a dictated amount of rows and use for loops for this case instead of a function code?

For your question, the polynomial functions given for temperautre (Ts1, Ts2, Tg1 and Tg2) are chosen based on the length (z). In the plots I have previously generated in excel, there's a significant trend difference onces z >=0.32, so I split the plot into two differenet trendlines, one for values z>=0 and z<=0.32 (where I want the code to chose functions Ts1 and Tg1) and for z>=0.32 and z<=0.41 (where I want to code to chose Ts2 and Tg2). In my mind this would generate a signle vector column that represent all the values for Ts and Tg along the desired length if that makes sense (basically replicate the excel plot). I'm trying to model a moving bed reactor, in which the other core equations (like density, reaction rate) depend on variation in temperature along the length of the reactor (z) so I'm trying to establish that part of the code first before moving on.

Thanks!

[Torsten]

What I meant is:

If you have a vector z with 2 components (in the above case z0 = [0 0.415]), do you want to evaluate your functions for both values of z independently and return a vector with two components for Tg and a vector of two components for Ts such that temp is 1x4 or 2x2 ?

[Jenna Broders]

Okay I understand now. I was thinking that that z "vector" would represents values spanning from 0 to 0.32 in this case. Similar to if you were creating a tspan using like t = linspace[0 0.32]. But in this scenario is just a 2 row vector of 0 and 0.32.

I'd ideally like the result to be something that would give two separate vecotrs, Ts and Tg, with them both being 1x2 vectors that were evaluated at the 2 values given in the z vector. So the first row would be the answer using the first number in the z vector and the second row would be the answer for the second number in th z vector

[Torsten]

z0 = [0 0.415];
[Ts,Tg] = feval(@temperature,z0);
Ts 
Tg

function [Ts,Tg] = temperature(z)
  Ts = zeros(size(z));
  Tg = zeros(size(z));
  idx1 = z>=0 & z <=0.32;
  idx2 = z>0.32;
  z1 = z(idx1);
  z2 = z(idx2);
  Ts1 = ((-3508440.83).*(z1.^5))+((3221859.34).*(z1.^4))-((1089760.82).*(z1.^3))+((165361.13).*(z1.^2))-((12015.96).*z1)+(1179.21);
  Ts2 = ((-1906981.2).*(z2.^3))+((2073664.15).*(z2.^2))-((754331.65).*(z2))+(92096.47);
  Tg1 = ((2796881.8).*(z1.^5))-((2553240.18).*(z1.^4))+((852324.51).*(z1.^3))-((124780.34).*(z1.^2))+((6092.92).*(z1))+(889.33);
  Tg2 = ((-1567143.11).*(z2.^3))+((1751141.7).*(z2.^2))-((653704.46).*(z2))+(81812.04);
  Ts(idx1) = Ts1;
  Ts(idx2) = Ts2;
  Tg(idx1) = Tg1;
  Tg(idx2) = Tg2;
end

[Jenna Broders]

Thanks for your help!