MATLAB: Partial Derivative and its roots

partial derivativeSymbolic Math Toolboxsyms

Dear Matlab Gurus

I want to calculate the symbolic partial derivative of a function w.r.t to 'a' and 'b'. The function is

    f=n*a*b*(1-b)*p*h^2*g^2)/(n*a*b*g^2*d*r^2+(1-a)*(1-b)*d*e*s^2)

Further more I also want to simplify and calculate the roots of both variable 'a' and 'b'. Actually I want to find the optimal point expression for both 'a' and 'b'. My code is

syms a
syms b
syms n
syms p
syms h
syms g
syms d
syms e
syms r
syms s
f=inline('(n*a*b*(1-b)*p*h^2*g^2)/(n*a*b*g^2*d*r^2+(1-a)*(1-b)*d*e*s^2)','a','b');
diff(diff(f(a,b),a),b)

however, the above code generates the following errors

Error using inlineeval (line 15)
Error in inline expression ==> (n*a*b*(1-b)*p*h^2*g^2)/(n*a*b*g^2*d*r^2+(1-a)*(1-b)*d*e*s^2)
Undefined function or variable 'n'.
Error in inline/subsref (line 24)
INLINE_OUT_ = inlineeval(INLINE_INPUTS_, INLINE_OBJ_.inputExpr, INLINE_OBJ_.expr);
 Error in test (line 12)
diff(diff(f(a,b),a),b)

Any kind of help is highly appreciated.

P.S: I have the symbolic math toolbox installed.

Regards

Best Answer

Hi Ali,

instead of the inline function you should use Anonymus Functions because inline will be removed in future releases.

syms a b n p h g d e r s
f=@(n,a,b,p,h,g,d,r,e,s) (n*a*b*(1-b)*p*h^2*g^2)/(n*a*b*g^2*d*r^2+(1-a)*(1-b)*d*e*s^2);
diff(diff(f,a),b)

Best regards,

Johannes

Related Solutions

MATLAB: Strange behavior of dsolve

This workaround works in R2019a:

clear all
close all
clc
syms g r s K P0 P(t)
assume(g==s)
eq1 = simplify(dsolve(diff(P,t) == r*P*(1-P/K)-s*P, P(0) == P0))
eq2 = simplify(dsolve(diff(P,t) == r*P*(1-P/K)-g*P, P(0) == P0))
res = isequal(eq1,eq2)

Setting the assumption below the dsolve command will fail:

clear all
close all
clc
syms g r s K P0 P(t)
eq1 = simplify(dsolve(diff(P,t) == r*P*(1-P/K)-s*P, P(0) == P0))
eq2 = simplify(dsolve(diff(P,t) == r*P*(1-P/K)-g*P, P(0) == P0))
assume(g==s)
res = isequal(eq1,eq2)

MATLAB: System of differential equations with many input arguments.

Try:

% l stands for \tau, u for \Theta, v for a, o for b (from utility) and j for ξ
syms S(t) P(t) c(t) h(t) z(t) A b d e g k l m n o r t u v w Y 
eqn1 = diff(S,t) == A*(b*S)^(1-g)*(z^g) - v*z - c - b*S*l;
eqn2 = diff(P,t) == u*z - d*P + j*c + (1-b)*S;
eqn3 = diff(h,t) ==  m*((b*S)-h);
eqn4 = diff(c,t) == (((b*(A*(1-g)*b^(1-g)*S^(-g)*z^g - l -r)) - ((A*g*(b*S)^(1-g)*z^(g-1) - v) * (1 - b + j*(d+r)) * (1/u))) * (u/(u - j*(v-g*A*(b*S)^(1-g)*z^(g-1)))))  +  ((1/(n*c^(-n)*(P^(-e)*b*S*h^(-o))^(1-n))) * ((exp(r*t)*w*m*b + b*(c*P^(-e)*h^(-o))^(1-n)*S^(-n)*b^(1-n)) + (j*e*P^(-e*(1-n)-1)*(c*b*S*h^(-o))^(1-n)) - ((1-n)*e*c^(-n)*P^(-e*(1-n)-1)*(b*S*h^(-o))^(1-n)*(u*z+(1-b)*S+j*c-d*P)) + (b*(1-n)*c^(-n)*(P^(-e)*h^(-o))^(1-n)*b^(1-n)*S^(-n)*(A*(b*S)^(1-g)*z^g - c - v*z - l*b*S))));
eqn5 = diff(z,t) == (((A*g*(b*S)^(1-g)*z^(g-1)-v) / (A*g*(g-1)*(b*S)^(1-g)*z^(g-2))) * ( d + ((e*P^(-e*(1-n)-1)*(c*b*S*h^(-o))^(1-n)*u) / (c^(-n)*(P*b*S*h^(-o))^(1-n))) - ((e*P^(-e*(1-n)-1)*(c*b*S*h^(-o))^(1-n)*j) / (c^(-n)*(P^(-e)*b*S*h^(-o))^(1-n))) + b*(A*(1-g)*b^(1-g)*S^(-g)*z^g - l + ((w*m*(u-j*(v-A*g*(b*S)^(1-g)*z^(g-1)))) / (exp(-r*t)*u*c^(-n)*(P^(-e)*b*S*h^(-o))^(1-n)))) + ((b^(1-n)*S^(-n)*(c*P^(-e)*h^(-o))^(1-n) * (u-j*(v-A*g*(b*S)^(1-g)*z^(g-1)))) / (u*e^(-n)*(P^(-e)*b*S*h^(-o))^(1-n))) + ((1-b)*(v-A*g*(b*S)^(1-g)*z^(g-1))*(1/u)))) - ((b*z/S) * (A*(b*S)^(1-g)*z^g - v*z - c - b*S*l));
[VF,Subs] = odeToVectorField(eqn1,eqn2,eqn3,eqn4,eqn5);
odefcn = matlabFunction(VF, 'Vars',{t,Y,A,b,d,e,g,k,l,m,n,o,r,u,v,w});
A=3;
b=0.6;
d=0.3;
e=0.4;
g=0.3;
k=4;
l=2;
m=0.7;
n=2;
o=0.6;
r=0.5;
u=10;	
v=2;
w=0.7;
tspan = [0 50];
y0 = [100 60 40 30 50];
[t,Y] = ode45(@(t,Y) odefcn(t,Y,A,b,d,e,g,k,l,m,n,o,r,u,v,w) , tspan, y0);
P = Y(:,1);	
S = Y(:,2);
h = Y(:,3);
c = Y(:,4);
z = Y(:,5);
x = b*S;		
q = A*((b*S).^(1-g)).*(z.^g);
figure(1)
plot(c,P)
title(sprintf('c-P'))
xlabel('Consumption (c)')
ylabel('Pollution (P)')
figure(2)
plot(c,x)
title(sprintf('C-x'))
xlabel('Consumption (c)')
ylabel('Recycling (x)')
figure(3)
plot(q,P)
title(sprintf('q-P'))
xlabel('Output (q)')
ylabel('Pollution (P)')
figure(4)
plot(q,x)
title(sprintf('q-x'))
xlabel('Output (q)')
ylabel('Recycling (x)')

3 reasons:

[t,Y] = ode45(@(t,Y) odefcn(t,A,Y,b,d,e,g,k,l,m,n,o,r,t,u,v,w) , tspan, y0);
%                             ^ ^                     ^
%                             | |                     |
%                             -----> Change order     ----> t already at first position
odefcn = matlabFunction(VF, 'Vars',{t,A,Y,b,d,e,g,k,l,m,n,o,r,u,v,w});
%                                   ^ ^                     
%                                   | |                     
%                  SAME HERE        -----> Change order

Best Answer

Related Solutions

MATLAB: Strange behavior of dsolve

MATLAB: System of differential equations with many input arguments.

Related Question