How do I pass out a vector of extra parameters using ODE15s

I'm writing up a code to solve a set of i differential equations and for each equation, I have three important parameters, which I need to export to the WS. I've tried to use evalin, since this was useful when I passed out a single variable.

Now I have to export a vector and I get the following error: "Error using evalin Undefined function 'aw_v' for input arguments of type 'double'."

This is part of the code:

for i=2:(par.nl-1)

          % Moisture content in the solids
          eqn2=@(phi_v)(1-par.epsilon)*rho_sv.*phi_v+par.epsilon*par.rho_a*((par.Psat/(0.0288*par.P))*(1-exp(-0.0004373*(par.T+57.926).*phi_v^(1.2626))))/(55.5556-((55.5556*par.Psat)/par.P)*(1-exp(-0.0004373*(par.T+57.926).*phi_v^(1.2626))))-W(i);

          phi_v(i)=fsolve(eqn2, par.phi); phi_v(i)=double(phi_v(i));

          % Water activity     
          aw_v(i)=isotherm(phi_v(i),par);

          % Humidity    
          H_v(i)=definition(aw_v(i),par);

          % Water vapour    
          WV(i)=par.epsilon*par.rho_a*H_v(i);

          dO(i)=(1/(par.epsilon.*par.Vt(i)))*(par.epsilon.*par.At(i)*par.D_O.*(O(i-1)-O(i))/par.dz-(1-par.epsilon)*par.rho_s.*par.Vt(i)*r_O-par.epsilon.*par.At(i+1)*par.D_O.*(O(i)-O(i+1))/par.dz);
          dC(i)=(1/(par.epsilon.*par.Vt(i)))*(par.epsilon.*par.At(i+1)*par.D_C*(C(i+1)-C(i))/par.dz+(1-par.epsilon)*par.rho_s.*par.Vt(i)*r_C-par.epsilon.*par.At(i)*par.D_C*(C(i)-C(i-1))/par.dz);
          dW(i)=(1./par.Vt(i))*(par.epsilon.*par.At(i+1)*par.D_W*(WV(i+1)-WV(i))*(1/par.epsilon)/par.dz-par.epsilon.*par.At(i)*par.D_W*(WV(i)-WV(i-1))*(1/par.epsilon)/par.dz);

          evalin('base', 'aw_v(i) = aw_v(i)'); evalin('base', 'H_v(i) = H_v(i)');  evalin('base', 'phi_v(i) = phi_v(i)');

      end

About this part: phi_v(i)=double(phi_v(i)); I find out that it is necessary since the outcome of the fsolve needs to be double for entering my function isotherm.

I really appreciate all your comments.

[Walter Roberson]

phi_v(i)=double(phi_v(i));

is not useful code.

Suppose somehow the data type of phi_v(1) were uint8 . Then you would have

phi_v(1) = double(phi_v(1))

The right hand side would create a double from the uint8, so you would have

phi_v(1) = some double

but the left hand side is datatype uint8 (in this scenario) so MATLAB will automatically convert the right hand side to match the left hand side, leading you to effectively

phi_v(1) = uint8( double(phi_v(1)) )

This is a useless operation. The only time it would have some potentially useful benefit is if phi_v were int64 or uint64 and for some reason you wanted to lose the extended precision in integers that int64 and uint64 offer compared to double.

If you want to change the data type of a vector or array you have to do the entire thing at one time, like

phi_v = double(phi_v);

In any case, the only possible output types of fsolve() are single or double.

Answer 1

[Walter Roberson]

assignin('base', 'aw_v(i)', aw_v(i));
assignin('base', 'H_v(i)', H_v(i));
assignin('base', 'phi_v(i)',  phi_v(i));

Beware: the results you will get in the base workspace will reflect whatever the final invocation of your ode function happened to calculate. The final invocation will not typically correspond to the final time. The ode functions do not evaluate only at the specified times: they evaluate at whatever times and boundary conditions they feel necessary in order to meet the integration tolerances, and they use the information to predict the results at the given times. If I remember correctly, ode15s will never move backwards in time, but some of the other routines like ode45 can move backwards in time.

If you need the results as of some particular time or boundary condition, then you are better off putting the calculation code into a function that is called by the ode function, and then calling the function "by hand" afterwards passing in the final time and boundary conditions output by the ode15s routine.

[Sara Jiménez]

Thank you for your quick reply and your advice!

I followed your first suggestion, but now I get the following error:

Error using assignin Invalid variable name "aw_v(i)" in ASSIGNIN.

[Walter Roberson]

You should probably remove the assignments from the loop and assign the entire vector after the loop.