How to Implement 3rd order equation in matlab script?
显示 更早的评论
I want to Implement this equation in matlab but I am not getting how to do it. Can anyone help me with this? psi1,psi2 and K1 are dependent on theta.
1 个评论
If both θ and Ψ are independent variables, then i is a function of two variables that produces a surface.
For example:
L = 160*membrane(1,100);
f = figure;
ax = axes;
s = surface(L); grid on
s.EdgeColor = 'none';
view(3)
ax.XLim = [1 201];
ax.YLim = [1 201];
ax.ZLim = [-53.4 160];
ax.CameraPosition = [-145.5 -229.7 283.6];
ax.CameraTarget = [77.4 60.2 63.9];
ax.CameraUpVector = [0 0 1];
ax.CameraViewAngle = 36.7;
l1 = light;
l1.Position = [160 400 80];
l1.Style = 'local';
l1.Color = [0 0.8 0.8];
l2 = light;
l2.Position = [.5 -1 .4];
l2.Color = [0.8 0.8 0];
s.FaceColor = [0.9 0.2 0.2];
xlabel('\theta', 'fontsize', 16), ylabel('\Psi', 'fontsize', 16), zlabel('i', 'fontsize', 16)
采纳的回答
Jon
2023-7-24
I would do it something like this:
function ival = ifun(theta,psi)
% Define constants
K2 = 25; % just for example, you put in actual value
K3 = 32.245; % just for example, you put in actual value
% Determine values of constants m and n based on values of psi
if psi > psi1(theta)
m = 1;
else
m = 0;
end
if psi > psi2(theta);
n = 1;
else
n = 0;
end
% Evaluate overall function
ival = K1(theta)*psi + m*K2*(psi - psi1(theta))^2 + n*K3*(psi - psi2(theta));
% Define helper functions
function k1val = K1(theta)
k1val = theta^2 + 3; % just for example you put in real function here
end
function psi1Val = psi1(theta)
psi1Val = 3*theta + theta^3; % just for example you put in real function here
end
function psi2Val = psi2(theta)
psi2Val = exp(theta) + sin(theta);% just for example you put in real function here
end
end
Then to actually evaluate ifun, either in a script or another function you would call it, for example
theta = pi/3;
psi = 0.9;
ival = ifun(theta,psi)
6 个评论
Prajwal Magadi
2023-7-25
Hello Jon,
Sorry for incomplete info, Actually psi1,psi2 and k1 are arrays dependent on theta. and K2=11,K3=185. after ivaluating I want to plot i vs psi.
Thank You
You could do it like this:
% Plot i as function of psi for fixed value of theta
theta = 20; % put your value of interest here
psi = linspace(0,1); % put your vector of psi values here
% loop to evaluate i at each value of psi
ivals = zeros(numel(psi),1); % preallocate
for k = 1:numel(psi)
ivals(k) = ifun(theta,psi(k));
end
% plot results
plot(psi,ivals)
xlabel('psi')
ylabel('i')
title(['i(theta,psi) for theta = ',num2str(theta)])
% Define function to evaluate i, could put in seperate m file, if you want to reuse
% it for other applications
function ival = ifun(theta,psi)
% Assign constants
K2 = 11; % just for example, you put in actual value
K3 = 185; % just for example, you put in actual value
% Assign theta dependent parameters (or perhaps read in from a data file if it can
% change)
thetaVals = 0:3:30;
K1Vals = [67,62.5,53.5,38,23.5,17,14,12,10,8.75,8];
psi1Vals = [0.25,0.25,0.25,0.175,0.2,0.225,0.335,0.46,0.47,0.485,0.485];
psi2Vals = [0.25,0.25,0.25,0.25,0.275,0.35,0.43,0.495,0.545,0.56,0.56];
% Look up values of theta dependent parameters
K1 = interp1(thetaVals,K1Vals,theta);
psi1 = interp1(thetaVals,psi1Vals,theta);
psi2 = interp1(thetaVals,psi2Vals,theta);
% Determine values of constants m and n based on values of psi
if psi > psi1
m = 1;
else
m = 0;
end
if psi > psi2
n = 1;
else
n = 0;
end
% Evaluate overall function
ival = K1*psi + m*K2*(psi - psi1)^2 + n*K3*(psi - psi2);
end
I see that you have theta in degrees. Usually physical equations are defined in terms of radians. You should check whether theta and psi should be in radians or degrees, and make sure that you are consistent.
Prajwal Magadi
2023-7-25
编辑：Prajwal Magadi
2023-7-25
Here is how you could solve for psi given i and plot results. If you don't actually have to solve for exact values of psi and just want the data plotted the other way you could just modigy the code I provided previously, swapping the arguments to plot and the xlabel and ylabel, to plot psi vs i instead of i vs psi.
% Plot psi as function of i for fixed value of theta
theta = 20; % put your value of interest here
ivals = 1:18; % range of i to be plotted
% loop to evaluate psi at each value of i
psi = zeros(numel(ivals),1); % preallocate
for k = 1:numel(ivals)
ival = ivals(k);
% Define function whose roots will give you values of psi
% need to do it inside of loop so that ival will be in scope
zfun = @(psi) ival - ifun(theta,psi);
psi(k) = fzero(zfun,[0,1]);
end
% plot results
figure
plot(ivals,psi)
ylabel('psi')
xlabel('i')
title(['psi(theta,i) for theta = ',num2str(theta)])
% Define function to evaluate i, could put in separate m file, if you want to reuse
% it for other applications
function ival = ifun(theta,psi)
% Assign constants
K2 = 11; % just for example, you put in actual value
K3 = 185; % just for example, you put in actual value
% Assign theta dependent parameters (or perhaps read in from a data file if it can
% change)
thetaVals = 0:3:30;
K1Vals = [67,62.5,53.5,38,23.5,17,14,12,10,8.75,8];
psi1Vals = [0.25,0.25,0.25,0.175,0.2,0.225,0.335,0.46,0.47,0.485,0.485];
psi2Vals = [0.25,0.25,0.25,0.25,0.275,0.35,0.43,0.495,0.545,0.56,0.56];
% Look up values of theta dependent parameters
K1 = interp1(thetaVals,K1Vals,theta);
psi1 = interp1(thetaVals,psi1Vals,theta);
psi2 = interp1(thetaVals,psi2Vals,theta);
% Determine values of constants m and n based on values of psi
if psi > psi1
m = 1;
else
m = 0;
end
if psi > psi2
n = 1;
else
n = 0;
end
% Evaluate overall function
ival = K1*psi + m*K2*(psi - psi1)^2 + n*K3*(psi - psi2);
end
Jon
2023-7-27
Did this answer your question? If so please accept the answer so others will know a solution is available. Otherwise what remaining issues do you have?
更多回答（1 个）
If 
and 
are functions of θ, and θ is the independent variable, then what exactly is Ψ? How does it look like on a graph? Can you sketch it?
and are functions of θ, and θ is the independent variable, then what exactly is Ψ? How does it look like on a graph? Can you sketch it?The i equation is not an issue, but it depends on m and n, and they are merely binary signals, consisting of only 1's and 0's arranged in some order. Technically, they are just Heaviside step functions that switch when the conditions are met.
See example for m:
theta = linspace(-1, 1, 2001);
Psi = theta;
Psi1 = sin(2*pi*Psi);
plot(theta, Psi, theta, Psi1), grid on
xlabel('\theta', 'fontsize', 16)
legend('\Psi', '\Psi_{1}', 'fontsize', 16, 'location', 'east')
m = 1/2*(sign(Psi - Psi1) + 1); % the math of Heaviside step function
plot(theta, m, 'linewidth', 2), grid on, ylim([-0.5 1.5])
xlabel('\theta', 'fontsize', 16)
ylabel('m', 'fontsize', 16)
I think you get the idea for now. But the real Ψ does not necessarily behave like a straight line.
类别
在 帮助中心 和 File Exchange 中查找有关 Characters and Strings 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
