필터 지우기
필터 지우기

how to use str2func inside functions that are solved by ODE45

조회 수: 2 (최근 30일)
Mehdi Gh
Mehdi Gh 2017년 3월 23일
댓글: Thomas Dixon 2020년 12월 3일
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.

채택된 답변

Jan
Jan 2017년 3월 23일
편집: Jan 2017년 3월 23일
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
Jan
Jan 2017년 3월 23일
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
Thomas Dixon 2020년 12월 3일
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

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Loops and Conditional Statements에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by