how to use str2func inside functions that are solved by ODE45
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
Mehdi Gh
am 23 Mär. 2017
Kommentiert: Thomas Dixon
am 3 Dez. 2020
I need to input functions as strings to the odefun that will be solved by ODE45. Problem is that the odefun takes the form of a cell and therefore can't be passed to ODE45 anymore. Is there any workaround this?
As an example:
function dc = pinene_kin2(t,c,k)
dc = cell(5,1); %5 components
str1 = '- k(1)*c(1) + k(2)*c(2)';
str2 = '+ k(1)*c(1) - k(2)*c(2)';
tmp1 = str2func(['@(c,k) ' str1]);
tmp2 = str2func(['@(c,k) ' str2]);
dc{1} = @(c)tmp1(c,k);
dc{2} = @(c)tmp2(c,k);
This function when called through ODE45:
ourparams = [.59 .3];
[T,C] = ode15s(@(t,c) pinene_kin2(t,c,our_params),[0 4],[100 0]);
produces the following error: Error using odearguments (line 110) Inputs must be floats, namely single or double.
NOTE: writing the equation as it is recommended by matlab doc is not acceptable. the odefun should be generated by the strings and then turned to the form required by ODE solver.
0 Kommentare
Akzeptierte Antwort
Jan
am 23 Mär. 2017
Bearbeitet: Jan
am 23 Mär. 2017
It is really weird to provide the functions as string. Most of all converting them to a function in each iteration is cruel.
If there is the need to create the equation as a string, why not creating an M-file dynamically, which contains this string and must be parsed once only for each integration?
Nevertheless, the problem of your code is something else: The function to be integrated must return a numerical array, not a cell:
function dc = pinene_kin2(t,c,k)
dc = zeros(2, 1);
str1 = '- k(1)*c(1) + k(2)*c(2)';
str2 = '+ k(1)*c(1) - k(2)*c(2)';
tmp1 = str2func(['@(c,k) ' str1]);
tmp2 = str2func(['@(c,k) ' str2]);
dc(1) = tmp1(c,k);
dc(2) = tmp2(c,k);
Why does your original function reply 5 components, although the initial values has only 2? Why did it reply a cell array containing anonymous functions? ODE15S requires numerical data. Do you want a symbolical integration?
Again: I cannot imagine that there can be a convincing reason to provide the functions as strings. It hurts to write the code I've suggest. Much better:
function dc = pinene_kin2(t,c,k)
dc = [- k(1)*c(1) + k(2)*c(2); ...
k(1)*c(1) - k(2)*c(2)];
and this M-file is a string also, when you write it to the hard disk.
3 Kommentare
Jan
am 23 Mär. 2017
Instead of providing the strings as inputs, which are converted to functions in each iteration, you can perform the conversion once and use the function handles as addition parameters.
Perhaps I've worked too long with students, who tried if my program can integrate functions like:
'system(''format C:'')'
(Of course this example does not work - unfortunately the students have been much smarter. :-) )
Thomas Dixon
am 3 Dez. 2020
Hi Jan,
I am doing something very similar and the option to have this write dynamic m-files to hard disk may/may not cause me a problem in my computation (HPC parallel process, I am not an expert but rule of thumb was limit read/write operations or hard disk space)
I am trying to do:
str2func(['@(x,A) ' dA])
where dA is:
+k(3)*A(3)*k(2)*conj(A(2))*exp(1i*(k(3)-k(2))*x) +k(4)*A(4)*k(3)*conj(A(3))*exp(1i*(k(4)-k(3))*x) +k(6)*A(6)*k(5)*conj(A(5))*exp(1i*(k(6)-k(5))*x); +k(3)*A(3)*k(1)*conj(A(1))*exp(1i*(k(3)-k(1))*x) +k(5)*A(5)*k(3)*conj(A(3))*exp(1i*(k(5)-k(3))*x) +k(6)*A(6)*k(4)*conj(A(4))*exp(1i*(k(6)-k(4))*x); -k(1)*A(1)*k(2)*A(2)*exp(1i*(k(1)+k(2))*x) +k(4)*A(4)*k(1)*conj(A(1))*exp(1i*(k(4)-k(1))*x) +k(5)*A(5)*k(2)*conj(A(2))*exp(1i*(k(5)-k(2))*x) +k(6)*A(6)*k(3)*conj(A(3))*exp(1i*(k(6)-k(3))*x); -k(1)*A(1)*k(3)*A(3)*exp(1i*(k(1)+k(3))*x) +k(6)*A(6)*k(2)*conj(A(2))*exp(1i*(k(6)-k(2))*x); -k(2)*A(2)*k(3)*A(3)*exp(1i*(k(2)+k(3))*x) +k(6)*A(6)*k(1)*conj(A(1))*exp(1i*(k(6)-k(1))*x); -k(1)*A(1)*k(5)*A(5)*exp(1i*(k(1)+k(5))*x) -k(2)*A(2)*k(4)*A(4)*exp(1i*(k(2)+k(4))*x) -0.5*k(3)*A(3)*k(3)*A(3)*exp(1i*(k(3)+k(3))*x)
I know I am not quite there yet with the basics. The manually constructed code is like this:
% dA = @(x,A)[-(maxBeta/2)*ks*ki*A(2)*A(3)*exp(1i*(ks+ki-kp)*x)+(maxBeta/2)*kshg*A(4)*kp*conj(A(1))*exp(1i*(kshg-2*kp)*x)+(maxBeta/2)*kps*ks*A(5)*conj(A(2))*exp(1i*(kps-ks-kp)*x)+(maxBeta/2)*kpi*ki*A(6)*conj(A(3))*exp(1i*(kpi-ki-kp)*x);
% (maxBeta/2)*ki*kp*conj(A(3))*A(1)*exp(1i*(kp-ks-ki)*x)+(maxBeta/2)*kp*kps*conj(A(1))*A(5)*exp(1i*(kps-kp-ks)*x)+(maxBeta/2)*kshg*kpi*conj(A(6))*A(4)*exp(1i*(kshg-kpi-ks)*x);
% (maxBeta/2)*ks*kp*conj(A(2))*A(1)*exp(1i*(kp-ks-ki)*x)+(maxBeta/2)*kp*kpi*conj(A(1))*A(6)*exp(1i*(kpi-kp-ki)*x)+(maxBeta/2)*kshg*kps*conj(A(5))*A(4)*exp(1i*(kshg-kps-ki)*x);
% -(maxBeta/4)*kp^2*A(1)^2*exp(1i*(2*kp-kshg)*x)-(maxBeta/2)*ki*kps*A(3)*A(5)*exp(1i*(ki+kps-kshg)*x)-(maxBeta/2)*kpi*ks*A(6)*A(2)*exp(1i*(kpi+ks-kshg)*x);
% -(maxBeta/2)*kp*ks*A(1)*A(2)*exp(1i*(kp+ks-kps)*x)+(maxBeta/2)*ki*kshg*conj(A(3))*A(4)*exp(1i*(kshg-ki-kps)*x);
% -(maxBeta/2)*kp*ki*A(1)*A(3)*exp(1i*(kp+ki-kpi)*x)+(maxBeta/2)*ks*kshg*conj(A(2))*A(4)*exp(1i*(kshg-ks-kpi)*x)];
%
which when given to:
As0=0.1;
Ap0=1;
[x,A] = ode45(dA, x ,[0 As0 Ap0 0 0 0]);
works well.
Any and all help is appreciated
Weitere Antworten (0)
Siehe auch
Kategorien
Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!