Where can I find MATLAB license files? Where can I find MATLAB license files?  Where can I find MATLAB license files?   MATLAB Answers — New Questions

​

Hi there,
I’m trying to create a matlab app in which users can create geometric objects and use them for some calculations later.
These are displayed in a 3D axis, rotated slightly. I am trying to create user functions to draw lines inside that axis. For this, I wan to implement continuously displayed x,y coordinates, or crosshairs or similar, as well as data snapping. Users will only draw 2D objects on the x,y plane. I still want users to be able to rotate the plot, zoom and pan normally with the mouse.
If I use a WindowButtonMotionFcn, then I lose plot interactivity with the mouse. If I use other means of getting the cursor position, like java.awt.MouseInfo.getPointerInfo().getLocation() and getpixelposition() or similar gui coordinate functions, I have to deal with a rotated camera and axis transformation. This makes "hit detection" on the data very cumbersome.
Using axes.CurrentPoint seems useless, because it only updates when the mouse is clicked, not continuously.
Any ideas on how to do this?Hi there,
I’m trying to create a matlab app in which users can create geometric objects and use them for some calculations later.
These are displayed in a 3D axis, rotated slightly. I am trying to create user functions to draw lines inside that axis. For this, I wan to implement continuously displayed x,y coordinates, or crosshairs or similar, as well as data snapping. Users will only draw 2D objects on the x,y plane. I still want users to be able to rotate the plot, zoom and pan normally with the mouse.
If I use a WindowButtonMotionFcn, then I lose plot interactivity with the mouse. If I use other means of getting the cursor position, like java.awt.MouseInfo.getPointerInfo().getLocation() and getpixelposition() or similar gui coordinate functions, I have to deal with a rotated camera and axis transformation. This makes "hit detection" on the data very cumbersome.
Using axes.CurrentPoint seems useless, because it only updates when the mouse is clicked, not continuously.
Any ideas on how to do this? Hi there,
I’m trying to create a matlab app in which users can create geometric objects and use them for some calculations later.
These are displayed in a 3D axis, rotated slightly. I am trying to create user functions to draw lines inside that axis. For this, I wan to implement continuously displayed x,y coordinates, or crosshairs or similar, as well as data snapping. Users will only draw 2D objects on the x,y plane. I still want users to be able to rotate the plot, zoom and pan normally with the mouse.
If I use a WindowButtonMotionFcn, then I lose plot interactivity with the mouse. If I use other means of getting the cursor position, like java.awt.MouseInfo.getPointerInfo().getLocation() and getpixelposition() or similar gui coordinate functions, I have to deal with a rotated camera and axis transformation. This makes "hit detection" on the data very cumbersome.
Using axes.CurrentPoint seems useless, because it only updates when the mouse is clicked, not continuously.
Any ideas on how to do this? appdesigner, app designer, windowbuttonmotionfcn, axes, matlab gui MATLAB Answers — New Questions

​

C-function code generation appears to work differently when generating code inside of a system composer software architecture model isntead of a standalone simulink model. I am trying to insert certain processer primitaves (mutexes) into my simulink model such that it is incorporated into the autogenerated code. It is a few simple lines of c++ code, so I opted to use the C-function block. In both cases, I have specified in the c-function dialogue window "Generate code as-is (optimizations off)"

Inside of an export-function simulink model, I created a c-function block that includes a comment and then a simple line of code, which generates correctly:

However, if I copy this exact same c-function into an export-function model within a system composer software architecture model, it no longer generates my my specified code. Instead, it generates function calls but does not define the functions anywhere. (Even doing a simple text search of my entire workspace only shows that the function call shown below gets generated).C-function code generation appears to work differently when generating code inside of a system composer software architecture model isntead of a standalone simulink model. I am trying to insert certain processer primitaves (mutexes) into my simulink model such that it is incorporated into the autogenerated code. It is a few simple lines of c++ code, so I opted to use the C-function block. In both cases, I have specified in the c-function dialogue window "Generate code as-is (optimizations off)"

Inside of an export-function simulink model, I created a c-function block that includes a comment and then a simple line of code, which generates correctly:

However, if I copy this exact same c-function into an export-function model within a system composer software architecture model, it no longer generates my my specified code. Instead, it generates function calls but does not define the functions anywhere. (Even doing a simple text search of my entire workspace only shows that the function call shown below gets generated). C-function code generation appears to work differently when generating code inside of a system composer software architecture model isntead of a standalone simulink model. I am trying to insert certain processer primitaves (mutexes) into my simulink model such that it is incorporated into the autogenerated code. It is a few simple lines of c++ code, so I opted to use the C-function block. In both cases, I have specified in the c-function dialogue window "Generate code as-is (optimizations off)"

Inside of an export-function simulink model, I created a c-function block that includes a comment and then a simple line of code, which generates correctly:

However, if I copy this exact same c-function into an export-function model within a system composer software architecture model, it no longer generates my my specified code. Instead, it generates function calls but does not define the functions anywhere. (Even doing a simple text search of my entire workspace only shows that the function call shown below gets generated). software architecture models, c-functions, export-function models, embedded coder, system composer MATLAB Answers — New Questions

​

The code below solves a 1-D, transient heat transfer problem set up as in general PDE format. The solution is plotted in color across the domain from 0 to 0.1 after 10 seconds have elapsed. What is the best way to plot the temperature across the length of this domain at this final time?

Thanks

clear all;

%% Create transient thermal model

thermalmodel = createpde(1);

R1= [3,4,0,0.1,0.1,0,0,0,1,1]’;
gd= [R1];

sf= ‘R1’;
ns = char(‘R1’);
ns = ns’;
dl = decsg(gd,sf,ns);

%% Create & plot geometry

geometryFromEdges(thermalmodel,dl);
pdegplot(thermalmodel,"EdgeLabels","on","FaceLabels","on")
xlim([0 0.1])
ylim([-1 1])
% axis equal

%% Generate and plot mesh

generateMesh(thermalmodel)
figure
pdemesh(thermalmodel)
title("Mesh with Quadratic Triangular Elements")

%% Apply BCs

% Edge 4 is left edge; Edge 2 is right side

applyBoundaryCondition(thermalmodel, "dirichlet",Edge=[4],u=100);
applyBoundaryCondition(thermalmodel, "dirichlet",Edge=[2],u=20);

%% Apply thermal properties [copper]

rho= 8933 %
cp= 385 %
rhocp= rho*cp %
k= 401 % W/mK

%% Define uniform volumetric heat generation rate

Qgen= 0 % W/m3

%% Define coefficients for generic Governing Equation to be solved

m= 0
a= 0
d= rhocp
c= [k]
f= [Qgen]

specifyCoefficients(thermalmodel, m=0, d=rhocp, c=k, a=0, f=f);

%% Apply initial condition

setInitialConditions(thermalmodel, 20);

%% Define time limits

tlist= 0: 1: 10;

thermalresults= solvepde(thermalmodel, tlist);

% Plot results

sol = thermalresults.NodalSolution;

subplot(2,2,1)
pdeplot(thermalmodel,"XYData",sol(:,11), …
"Contour","on",…
"ColorMap","jet")The code below solves a 1-D, transient heat transfer problem set up as in general PDE format. The solution is plotted in color across the domain from 0 to 0.1 after 10 seconds have elapsed. What is the best way to plot the temperature across the length of this domain at this final time?

Thanks

clear all;

%% Create transient thermal model

thermalmodel = createpde(1);

R1= [3,4,0,0.1,0.1,0,0,0,1,1]’;
gd= [R1];

sf= ‘R1’;
ns = char(‘R1’);
ns = ns’;
dl = decsg(gd,sf,ns);

%% Create & plot geometry

geometryFromEdges(thermalmodel,dl);
pdegplot(thermalmodel,"EdgeLabels","on","FaceLabels","on")
xlim([0 0.1])
ylim([-1 1])
% axis equal

%% Generate and plot mesh

generateMesh(thermalmodel)
figure
pdemesh(thermalmodel)
title("Mesh with Quadratic Triangular Elements")

%% Apply BCs

% Edge 4 is left edge; Edge 2 is right side

applyBoundaryCondition(thermalmodel, "dirichlet",Edge=[4],u=100);
applyBoundaryCondition(thermalmodel, "dirichlet",Edge=[2],u=20);

%% Apply thermal properties [copper]

rho= 8933 %
cp= 385 %
rhocp= rho*cp %
k= 401 % W/mK

%% Define uniform volumetric heat generation rate

Qgen= 0 % W/m3

%% Define coefficients for generic Governing Equation to be solved

m= 0
a= 0
d= rhocp
c= [k]
f= [Qgen]

specifyCoefficients(thermalmodel, m=0, d=rhocp, c=k, a=0, f=f);

%% Apply initial condition

setInitialConditions(thermalmodel, 20);

%% Define time limits

tlist= 0: 1: 10;

thermalresults= solvepde(thermalmodel, tlist);

% Plot results

sol = thermalresults.NodalSolution;

subplot(2,2,1)
pdeplot(thermalmodel,"XYData",sol(:,11), …
"Contour","on",…
"ColorMap","jet") The code below solves a 1-D, transient heat transfer problem set up as in general PDE format. The solution is plotted in color across the domain from 0 to 0.1 after 10 seconds have elapsed. What is the best way to plot the temperature across the length of this domain at this final time?

Thanks

clear all;

%% Create transient thermal model

thermalmodel = createpde(1);

R1= [3,4,0,0.1,0.1,0,0,0,1,1]’;
gd= [R1];

sf= ‘R1’;
ns = char(‘R1’);
ns = ns’;
dl = decsg(gd,sf,ns);

%% Create & plot geometry

geometryFromEdges(thermalmodel,dl);
pdegplot(thermalmodel,"EdgeLabels","on","FaceLabels","on")
xlim([0 0.1])
ylim([-1 1])
% axis equal

%% Generate and plot mesh

generateMesh(thermalmodel)
figure
pdemesh(thermalmodel)
title("Mesh with Quadratic Triangular Elements")

%% Apply BCs

% Edge 4 is left edge; Edge 2 is right side

applyBoundaryCondition(thermalmodel, "dirichlet",Edge=[4],u=100);
applyBoundaryCondition(thermalmodel, "dirichlet",Edge=[2],u=20);

%% Apply thermal properties [copper]

rho= 8933 %
cp= 385 %
rhocp= rho*cp %
k= 401 % W/mK

%% Define uniform volumetric heat generation rate

Qgen= 0 % W/m3

%% Define coefficients for generic Governing Equation to be solved

m= 0
a= 0
d= rhocp
c= [k]
f= [Qgen]

specifyCoefficients(thermalmodel, m=0, d=rhocp, c=k, a=0, f=f);

%% Apply initial condition

setInitialConditions(thermalmodel, 20);

%% Define time limits

tlist= 0: 1: 10;

thermalresults= solvepde(thermalmodel, tlist);

% Plot results

sol = thermalresults.NodalSolution;

subplot(2,2,1)
pdeplot(thermalmodel,"XYData",sol(:,11), …
"Contour","on",…
"ColorMap","jet") plotting MATLAB Answers — New Questions

​

Hello,
Is the code generation tunable with a transition variant, set with "Variant Activation Time = Code compile" with r2024b?
Indeed, I have a compilation warning with my generated code. The code generation creates a static function with #if inside the function:

This way, the compilation gives a compilation warning saying that there is a function defined but not used.
If you have a solution, I am interested.
Regards,
JeanHello,
Is the code generation tunable with a transition variant, set with "Variant Activation Time = Code compile" with r2024b?
Indeed, I have a compilation warning with my generated code. The code generation creates a static function with #if inside the function:

This way, the compilation gives a compilation warning saying that there is a function defined but not used.
If you have a solution, I am interested.
Regards,
Jean Hello,
Is the code generation tunable with a transition variant, set with "Variant Activation Time = Code compile" with r2024b?
Indeed, I have a compilation warning with my generated code. The code generation creates a static function with #if inside the function:

This way, the compilation gives a compilation warning saying that there is a function defined but not used.
If you have a solution, I am interested.
Regards,
Jean variant, code generation MATLAB Answers — New Questions

​

The model:

<</matlabcentral/answers/uploaded_files/27438/3.PNG>>

The simMechanics system is as follow:

<</matlabcentral/answers/uploaded_files/27436/1.PNG>>

simin is input from workspace and is the motion of Revolute Joint as the time goes on.

<</matlabcentral/answers/uploaded_files/27437/2.PNG>>

But it shows some kind of error:
In the dynamically coupled component containing Revolute Joint Revolute_Joint, there are fewer joint primitive degrees of freedom with automatically computed force or torque (0) than with motion from inputs (1). The prescribed motion trajectories in this component may not be achievable. Solve this problem by reducing the number of joint primitives with motion from inputs or increasing the number of joint primitives with automatically computed force or torque. Resolve this issue in order to simulate the model.

Does anybody knows what’s wrong and what should I do to solve the problem.The model:

<</matlabcentral/answers/uploaded_files/27438/3.PNG>>

The simMechanics system is as follow:

<</matlabcentral/answers/uploaded_files/27436/1.PNG>>

simin is input from workspace and is the motion of Revolute Joint as the time goes on.

<</matlabcentral/answers/uploaded_files/27437/2.PNG>>

But it shows some kind of error:
In the dynamically coupled component containing Revolute Joint Revolute_Joint, there are fewer joint primitive degrees of freedom with automatically computed force or torque (0) than with motion from inputs (1). The prescribed motion trajectories in this component may not be achievable. Solve this problem by reducing the number of joint primitives with motion from inputs or increasing the number of joint primitives with automatically computed force or torque. Resolve this issue in order to simulate the model.

Does anybody knows what’s wrong and what should I do to solve the problem. The model:

<</matlabcentral/answers/uploaded_files/27438/3.PNG>>

The simMechanics system is as follow:

<</matlabcentral/answers/uploaded_files/27436/1.PNG>>

simin is input from workspace and is the motion of Revolute Joint as the time goes on.

<</matlabcentral/answers/uploaded_files/27437/2.PNG>>

But it shows some kind of error:
In the dynamically coupled component containing Revolute Joint Revolute_Joint, there are fewer joint primitive degrees of freedom with automatically computed force or torque (0) than with motion from inputs (1). The prescribed motion trajectories in this component may not be achievable. Solve this problem by reducing the number of joint primitives with motion from inputs or increasing the number of joint primitives with automatically computed force or torque. Resolve this issue in order to simulate the model.

Does anybody knows what’s wrong and what should I do to solve the problem. simulink, simmechanics MATLAB Answers — New Questions

​

I encountered a nonlinear fitting problem, and the fitting formula is shown in Equation (1), which includes two infinite integrals (in practice, the integration range can be set from 0.01E-6 to 200E-6).

In these formulas, except for x and y being vectors, all other variables are scalars, and Rmedian and sigma are the parameters to be fitted.

I found a related post and tried to write the code based on it, but it keeps reporting errors. The error message seems to indicate that the vector dimensions are inconsistent, preventing the operation from proceeding. However, these functions are all calculations for individual scalars.
Error using /
Matrix dimensions must agree.

Error in dsdmain>@(r)1/(2*r*sigma*sqrt(2*pi))*exp(-(log(2*r)-log(2*Rmean)^2)/(2*sigma^2)) (line 13)
gauss = @(r) 1/(2*r*sigma*sqrt(2*pi))* exp( -(log(2*r)-log(2*Rmean)^2)/(2*sigma^2) );
My question is: Can I refer to the content of this post to solve my problem? If yes, what does this error message mean? If not, how should I resolve my problem? (Note: The range of Rmedian is 1E-6 to 5E-6)

After modifying the code according to Walter Roberson and Torsten’s suggestion, the program no longer throws errors. But no matter how I set the initial values on my end, it always prompts:
Initial point is a local minimum.
Optimization completed because the size of the gradient at the initial point
is less than the default value of the optimality tolerance.
<stopping criteria details>
theta = initial point
1.0e-05 *
0.2000
0.1000
Optimization completed: The final point is the initial point.
The first-order optimality measure, 0.000000e+00, is less than
options.OptimalityTolerance = 1.000000e-06.
Optimization Metric Options
relative first-order optimality = 0.00e+00 OptimalityTolerance = 1e-06 (default)

I have checked all the formulas and the units of the variables, and I didn’t find any problems.

————————————– Beblow is my code for issue reproduction ———————————————————-
function testmain()
clc
function Kvec = model(param,xdata)
% Vector of Kd for every xdata:
Kvec = zeros(size(xdata));
Rmean = param(1);
Rstd = param(2);

for i = 1:length(xdata)
model = @(r) unified(xdata(i),r,Rmean,Rstd,delta,Delta,D,lambda);
Kvec(i) = integral(model,0.01E-6,200E-6);
end
end
function s = unified(g,R,Rmean,Rstd,delta,Delta,D,lambda)
%unified Unified fitting model for DSD
% exponentional combined

factor = 1./(2*R*Rstd*sqrt(2*pi)) *2 ; % int(P(r)) = 0.5,1/0.5=2
p1 = -(log(2*R)-log(2*Rmean)).^2/(2*Rstd^2);
c = -2*2.675E8^2*g.^2/D;

tmp = 0;
for il = 1:length(lambda)
a2 = (lambda(il)./R).^2;
an4 = (R/lambda(il)).^4;
Psi = 2+exp(-a2*D*(Delta-delta))-2*exp(-a2*D*delta)-2*exp(-a2*D*Delta)+exp(-a2*D*(Delta+delta));
tmp = tmp+an4./(a2.*R.*R-2).*(2*delta-Psi./(a2*D));
end
p2 = c*tmp;

s = factor.*exp(p1+p2);
end
Delta = 0.075;
delta = 0.002;
D = 0.098E-9;
lambda = [2.0816 5.9404 9.2058 12.4044 15.5792 18.7426 21.8997 25.0528 28.2034 31.3521];

g = [ 0.300616, 0.53884, 0.771392, 1.009616, 1.24784, 1.480392, 1.718616, 1.95684, 2.189392, 2.427616, 2.66584, 2.898392 ];
xdata = g;
ydata = [100, 91.16805426, 80.52955192, 67.97705378, 55.1009735,41.87307917, 30.39638776, 21.13515607, 13.7125649, 8.33083767, 5.146756077, 2.79768609];
ydata = ydata/ydata(1); % normalize

% Inital guess for parameters:
Rmean0 = 2E-6;
Rstd0 = 1E-6;
p0 = [Rmean0;Rstd0];
% lsqcurvefit is in the optimization toolbox.
% fit, from the curve fitting toolbox may be an alternative
theta = lsqcurvefit(@model,p0,xdata,ydata,[0.01E-6;0.1E-6],[10E-6,2E-6])
endI encountered a nonlinear fitting problem, and the fitting formula is shown in Equation (1), which includes two infinite integrals (in practice, the integration range can be set from 0.01E-6 to 200E-6).

In these formulas, except for x and y being vectors, all other variables are scalars, and Rmedian and sigma are the parameters to be fitted.

I found a related post and tried to write the code based on it, but it keeps reporting errors. The error message seems to indicate that the vector dimensions are inconsistent, preventing the operation from proceeding. However, these functions are all calculations for individual scalars.
Error using /
Matrix dimensions must agree.

Error in dsdmain>@(r)1/(2*r*sigma*sqrt(2*pi))*exp(-(log(2*r)-log(2*Rmean)^2)/(2*sigma^2)) (line 13)
gauss = @(r) 1/(2*r*sigma*sqrt(2*pi))* exp( -(log(2*r)-log(2*Rmean)^2)/(2*sigma^2) );
My question is: Can I refer to the content of this post to solve my problem? If yes, what does this error message mean? If not, how should I resolve my problem? (Note: The range of Rmedian is 1E-6 to 5E-6)

After modifying the code according to Walter Roberson and Torsten’s suggestion, the program no longer throws errors. But no matter how I set the initial values on my end, it always prompts:
Initial point is a local minimum.
Optimization completed because the size of the gradient at the initial point
is less than the default value of the optimality tolerance.
<stopping criteria details>
theta = initial point
1.0e-05 *
0.2000
0.1000
Optimization completed: The final point is the initial point.
The first-order optimality measure, 0.000000e+00, is less than
options.OptimalityTolerance = 1.000000e-06.
Optimization Metric Options
relative first-order optimality = 0.00e+00 OptimalityTolerance = 1e-06 (default)

I have checked all the formulas and the units of the variables, and I didn’t find any problems.

————————————– Beblow is my code for issue reproduction ———————————————————-
function testmain()
clc
function Kvec = model(param,xdata)
% Vector of Kd for every xdata:
Kvec = zeros(size(xdata));
Rmean = param(1);
Rstd = param(2);

for i = 1:length(xdata)
model = @(r) unified(xdata(i),r,Rmean,Rstd,delta,Delta,D,lambda);
Kvec(i) = integral(model,0.01E-6,200E-6);
end
end
function s = unified(g,R,Rmean,Rstd,delta,Delta,D,lambda)
%unified Unified fitting model for DSD
% exponentional combined

factor = 1./(2*R*Rstd*sqrt(2*pi)) *2 ; % int(P(r)) = 0.5,1/0.5=2
p1 = -(log(2*R)-log(2*Rmean)).^2/(2*Rstd^2);
c = -2*2.675E8^2*g.^2/D;

tmp = 0;
for il = 1:length(lambda)
a2 = (lambda(il)./R).^2;
an4 = (R/lambda(il)).^4;
Psi = 2+exp(-a2*D*(Delta-delta))-2*exp(-a2*D*delta)-2*exp(-a2*D*Delta)+exp(-a2*D*(Delta+delta));
tmp = tmp+an4./(a2.*R.*R-2).*(2*delta-Psi./(a2*D));
end
p2 = c*tmp;

s = factor.*exp(p1+p2);
end
Delta = 0.075;
delta = 0.002;
D = 0.098E-9;
lambda = [2.0816 5.9404 9.2058 12.4044 15.5792 18.7426 21.8997 25.0528 28.2034 31.3521];

g = [ 0.300616, 0.53884, 0.771392, 1.009616, 1.24784, 1.480392, 1.718616, 1.95684, 2.189392, 2.427616, 2.66584, 2.898392 ];
xdata = g;
ydata = [100, 91.16805426, 80.52955192, 67.97705378, 55.1009735,41.87307917, 30.39638776, 21.13515607, 13.7125649, 8.33083767, 5.146756077, 2.79768609];
ydata = ydata/ydata(1); % normalize

% Inital guess for parameters:
Rmean0 = 2E-6;
Rstd0 = 1E-6;
p0 = [Rmean0;Rstd0];
% lsqcurvefit is in the optimization toolbox.
% fit, from the curve fitting toolbox may be an alternative
theta = lsqcurvefit(@model,p0,xdata,ydata,[0.01E-6;0.1E-6],[10E-6,2E-6])
end I encountered a nonlinear fitting problem, and the fitting formula is shown in Equation (1), which includes two infinite integrals (in practice, the integration range can be set from 0.01E-6 to 200E-6).

In these formulas, except for x and y being vectors, all other variables are scalars, and Rmedian and sigma are the parameters to be fitted.

I found a related post and tried to write the code based on it, but it keeps reporting errors. The error message seems to indicate that the vector dimensions are inconsistent, preventing the operation from proceeding. However, these functions are all calculations for individual scalars.
Error using /
Matrix dimensions must agree.

Error in dsdmain>@(r)1/(2*r*sigma*sqrt(2*pi))*exp(-(log(2*r)-log(2*Rmean)^2)/(2*sigma^2)) (line 13)
gauss = @(r) 1/(2*r*sigma*sqrt(2*pi))* exp( -(log(2*r)-log(2*Rmean)^2)/(2*sigma^2) );
My question is: Can I refer to the content of this post to solve my problem? If yes, what does this error message mean? If not, how should I resolve my problem? (Note: The range of Rmedian is 1E-6 to 5E-6)

After modifying the code according to Walter Roberson and Torsten’s suggestion, the program no longer throws errors. But no matter how I set the initial values on my end, it always prompts:
Initial point is a local minimum.
Optimization completed because the size of the gradient at the initial point
is less than the default value of the optimality tolerance.
<stopping criteria details>
theta = initial point
1.0e-05 *
0.2000
0.1000
Optimization completed: The final point is the initial point.
The first-order optimality measure, 0.000000e+00, is less than
options.OptimalityTolerance = 1.000000e-06.
Optimization Metric Options
relative first-order optimality = 0.00e+00 OptimalityTolerance = 1e-06 (default)

I have checked all the formulas and the units of the variables, and I didn’t find any problems.

————————————– Beblow is my code for issue reproduction ———————————————————-
function testmain()
clc
function Kvec = model(param,xdata)
% Vector of Kd for every xdata:
Kvec = zeros(size(xdata));
Rmean = param(1);
Rstd = param(2);

for i = 1:length(xdata)
model = @(r) unified(xdata(i),r,Rmean,Rstd,delta,Delta,D,lambda);
Kvec(i) = integral(model,0.01E-6,200E-6);
end
end
function s = unified(g,R,Rmean,Rstd,delta,Delta,D,lambda)
%unified Unified fitting model for DSD
% exponentional combined

factor = 1./(2*R*Rstd*sqrt(2*pi)) *2 ; % int(P(r)) = 0.5,1/0.5=2
p1 = -(log(2*R)-log(2*Rmean)).^2/(2*Rstd^2);
c = -2*2.675E8^2*g.^2/D;

tmp = 0;
for il = 1:length(lambda)
a2 = (lambda(il)./R).^2;
an4 = (R/lambda(il)).^4;
Psi = 2+exp(-a2*D*(Delta-delta))-2*exp(-a2*D*delta)-2*exp(-a2*D*Delta)+exp(-a2*D*(Delta+delta));
tmp = tmp+an4./(a2.*R.*R-2).*(2*delta-Psi./(a2*D));
end
p2 = c*tmp;

s = factor.*exp(p1+p2);
end
Delta = 0.075;
delta = 0.002;
D = 0.098E-9;
lambda = [2.0816 5.9404 9.2058 12.4044 15.5792 18.7426 21.8997 25.0528 28.2034 31.3521];

g = [ 0.300616, 0.53884, 0.771392, 1.009616, 1.24784, 1.480392, 1.718616, 1.95684, 2.189392, 2.427616, 2.66584, 2.898392 ];
xdata = g;
ydata = [100, 91.16805426, 80.52955192, 67.97705378, 55.1009735,41.87307917, 30.39638776, 21.13515607, 13.7125649, 8.33083767, 5.146756077, 2.79768609];
ydata = ydata/ydata(1); % normalize

% Inital guess for parameters:
Rmean0 = 2E-6;
Rstd0 = 1E-6;
p0 = [Rmean0;Rstd0];
% lsqcurvefit is in the optimization toolbox.
% fit, from the curve fitting toolbox may be an alternative
theta = lsqcurvefit(@model,p0,xdata,ydata,[0.01E-6;0.1E-6],[10E-6,2E-6])
end curve fitting, integral MATLAB Answers — New Questions

​

Hello,
I have a question regarding the relationship between Simulink sample rate and HDL Coder target frequency.
Is it okay if the Simulink sample rate is equal to the target frequency specified in HDL Coder?
I’ve often heard that the hardware clock rate (target frequency) should generally be at least 2× the sample rate to ensure proper timing and pipelining. But in Simulink and HDL Coder, I’m unsure how strict this rule is.
Specifically:
If I set my Simulink sample rate to 120 MSPS and also specify 120 MHz as the HDL Coder target frequency, is this configuration valid?
Or should the target frequency always be greater than the sample rate, even in fully pipelined or parallel architectures?
I’d really appreciate any clarification on how these two values interact and whether matching them could lead to timing issues during synthesis or implementation.
Thanks in advance!Hello,
I have a question regarding the relationship between Simulink sample rate and HDL Coder target frequency.
Is it okay if the Simulink sample rate is equal to the target frequency specified in HDL Coder?
I’ve often heard that the hardware clock rate (target frequency) should generally be at least 2× the sample rate to ensure proper timing and pipelining. But in Simulink and HDL Coder, I’m unsure how strict this rule is.
Specifically:
If I set my Simulink sample rate to 120 MSPS and also specify 120 MHz as the HDL Coder target frequency, is this configuration valid?
Or should the target frequency always be greater than the sample rate, even in fully pipelined or parallel architectures?
I’d really appreciate any clarification on how these two values interact and whether matching them could lead to timing issues during synthesis or implementation.
Thanks in advance! Hello,
I have a question regarding the relationship between Simulink sample rate and HDL Coder target frequency.
Is it okay if the Simulink sample rate is equal to the target frequency specified in HDL Coder?
I’ve often heard that the hardware clock rate (target frequency) should generally be at least 2× the sample rate to ensure proper timing and pipelining. But in Simulink and HDL Coder, I’m unsure how strict this rule is.
Specifically:
If I set my Simulink sample rate to 120 MSPS and also specify 120 MHz as the HDL Coder target frequency, is this configuration valid?
Or should the target frequency always be greater than the sample rate, even in fully pipelined or parallel architectures?
I’d really appreciate any clarification on how these two values interact and whether matching them could lead to timing issues during synthesis or implementation.
Thanks in advance! simulink, sample rate, hdl coder, target frequency MATLAB Answers — New Questions

​

Hi everyone,
I’m trying to generate a uniform grid where the red lines intersect at a 45-degree angle to the horizontal lines in blue using the finite difference method in MATLAB.
However, the resulting grid is slightly off, as shown in the attached image. Intesection angle is uniform but higher than 45-degree.
I set computational domain as and a grid is generated on physical domain with coordinates .
As I want to generate a mesh described above, I aimed to solve the following PDE :

where we have (horizontal lines). The boundary conditions are .
The above equation comes down to the following:

Based on the above equation, I discretized to finite difference form (central difference for and forward differece for ).
Plus, I found the grid collapses when I set different angles (for example 90-degree (replace )).
Am I missing something in implementation?

I need your help.
Here is my code.
% Parameter settings
eta_max = 10;
xi_max = 10;
Delta_eta = 0.1;
Delta_xi = 0.1;
tol = 1e-6;
max_iterations = 100000;

% Grid size
eta_steps = floor(eta_max / Delta_eta) + 1;
xi_steps = floor(xi_max / Delta_xi) + 1;

% Initial conditions
x = zeros(eta_steps, xi_steps);
y = repmat((0:eta_steps-1)’ * Delta_eta, 1, xi_steps);

% Boundary conditions
x(1, 🙂 = linspace(0, xi_max, xi_steps);

% Finite difference method
converged = false;
iteration = 0;

while ~converged && iteration < max_iterations
iteration = iteration + 1;
max_error = 0;

for n = 1:eta_steps – 1
for m = 2:xi_steps – 1
x_xi = (x(n, m+1) – x(n, m-1)) / (2 * Delta_xi);
y_xi = (y(n, m+1) – y(n, m-1)) / (2 * Delta_xi);

new_x = x(n, m) + Delta_eta * ((pi/4) – y_xi)/(x_xi);
error_x = abs(new_x – x(n, m));
max_error = max([max_error, error_x]);

x(n + 1, m) = new_x;
end

x(n+1, xi_steps) = x(n+1, xi_steps-1) + Delta_xi;
x(n+1, 1) = x(n+1, 2) – Delta_xi;
end

if max_error < tol
converged = true;
end
end

% Plot the results
figure; hold on; axis equal; grid on;

for m = 1:xi_steps
plot(x(:, m), y(:, m), ‘r’);
end

for n = 1:eta_steps
plot(x(n, :), y(n, :), ‘b’);
end

title([‘Iterations until convergence: ‘, num2str(iteration)]);
xlabel(‘x’); ylabel(‘y’);Hi everyone,
I’m trying to generate a uniform grid where the red lines intersect at a 45-degree angle to the horizontal lines in blue using the finite difference method in MATLAB.
However, the resulting grid is slightly off, as shown in the attached image. Intesection angle is uniform but higher than 45-degree.
I set computational domain as and a grid is generated on physical domain with coordinates .
As I want to generate a mesh described above, I aimed to solve the following PDE :

where we have (horizontal lines). The boundary conditions are .
The above equation comes down to the following:

Based on the above equation, I discretized to finite difference form (central difference for and forward differece for ).
Plus, I found the grid collapses when I set different angles (for example 90-degree (replace )).
Am I missing something in implementation?

I need your help.
Here is my code.
% Parameter settings
eta_max = 10;
xi_max = 10;
Delta_eta = 0.1;
Delta_xi = 0.1;
tol = 1e-6;
max_iterations = 100000;

% Grid size
eta_steps = floor(eta_max / Delta_eta) + 1;
xi_steps = floor(xi_max / Delta_xi) + 1;

% Initial conditions
x = zeros(eta_steps, xi_steps);
y = repmat((0:eta_steps-1)’ * Delta_eta, 1, xi_steps);

% Boundary conditions
x(1, 🙂 = linspace(0, xi_max, xi_steps);

% Finite difference method
converged = false;
iteration = 0;

while ~converged && iteration < max_iterations
iteration = iteration + 1;
max_error = 0;

for n = 1:eta_steps – 1
for m = 2:xi_steps – 1
x_xi = (x(n, m+1) – x(n, m-1)) / (2 * Delta_xi);
y_xi = (y(n, m+1) – y(n, m-1)) / (2 * Delta_xi);

new_x = x(n, m) + Delta_eta * ((pi/4) – y_xi)/(x_xi);
error_x = abs(new_x – x(n, m));
max_error = max([max_error, error_x]);

x(n + 1, m) = new_x;
end

x(n+1, xi_steps) = x(n+1, xi_steps-1) + Delta_xi;
x(n+1, 1) = x(n+1, 2) – Delta_xi;
end

if max_error < tol
converged = true;
end
end

% Plot the results
figure; hold on; axis equal; grid on;

for m = 1:xi_steps
plot(x(:, m), y(:, m), ‘r’);
end

for n = 1:eta_steps
plot(x(n, :), y(n, :), ‘b’);
end

title([‘Iterations until convergence: ‘, num2str(iteration)]);
xlabel(‘x’); ylabel(‘y’); Hi everyone,
I’m trying to generate a uniform grid where the red lines intersect at a 45-degree angle to the horizontal lines in blue using the finite difference method in MATLAB.
However, the resulting grid is slightly off, as shown in the attached image. Intesection angle is uniform but higher than 45-degree.
I set computational domain as and a grid is generated on physical domain with coordinates .
As I want to generate a mesh described above, I aimed to solve the following PDE :

where we have (horizontal lines). The boundary conditions are .
The above equation comes down to the following:

Based on the above equation, I discretized to finite difference form (central difference for and forward differece for ).
Plus, I found the grid collapses when I set different angles (for example 90-degree (replace )).
Am I missing something in implementation?

I need your help.
Here is my code.
% Parameter settings
eta_max = 10;
xi_max = 10;
Delta_eta = 0.1;
Delta_xi = 0.1;
tol = 1e-6;
max_iterations = 100000;

% Grid size
eta_steps = floor(eta_max / Delta_eta) + 1;
xi_steps = floor(xi_max / Delta_xi) + 1;

% Initial conditions
x = zeros(eta_steps, xi_steps);
y = repmat((0:eta_steps-1)’ * Delta_eta, 1, xi_steps);

% Boundary conditions
x(1, 🙂 = linspace(0, xi_max, xi_steps);

% Finite difference method
converged = false;
iteration = 0;

while ~converged && iteration < max_iterations
iteration = iteration + 1;
max_error = 0;

for n = 1:eta_steps – 1
for m = 2:xi_steps – 1
x_xi = (x(n, m+1) – x(n, m-1)) / (2 * Delta_xi);
y_xi = (y(n, m+1) – y(n, m-1)) / (2 * Delta_xi);

new_x = x(n, m) + Delta_eta * ((pi/4) – y_xi)/(x_xi);
error_x = abs(new_x – x(n, m));
max_error = max([max_error, error_x]);

x(n + 1, m) = new_x;
end

x(n+1, xi_steps) = x(n+1, xi_steps-1) + Delta_xi;
x(n+1, 1) = x(n+1, 2) – Delta_xi;
end

if max_error < tol
converged = true;
end
end

% Plot the results
figure; hold on; axis equal; grid on;

for m = 1:xi_steps
plot(x(:, m), y(:, m), ‘r’);
end

for n = 1:eta_steps
plot(x(n, :), y(n, :), ‘b’);
end

title([‘Iterations until convergence: ‘, num2str(iteration)]);
xlabel(‘x’); ylabel(‘y’); finite difference method, grid generation, mesh generation, numerical analysis, differential equations MATLAB Answers — New Questions

​

Why numerical solutions is not consistent with the exact solutions? Also the error is too large. How to modify the code to get consistent results? See the attached image for the equation of motions and their exact solutions.
clear; clc;
tic;
c1 = 0; % real constant
c2 = 0; % real constant
a = 1; % real constant
b=0.1; % two values b=0.1, b=0.01
lambda = -1./(2*(1 + 1i)); % complex or imaginary
N = 1000;
X = linspace(-10, 10, N);
T = linspace(0, 10, N);

% Calculate initial conditions
[x, ~] = meshgrid(X, 0); % t = 0 for initial conditions
%psi_init = A * exp(-1i * lambda * x);
%phi_init = B * exp(1i * lambda * x);
nq_init = -a*(b.^5*(1+(0 + c2).^2).^2 – 8*a*(0+c2).*(x+c1) – 4*a*b.^4*(0+c2).*(1+(0+c2).^2).*(x+c1) – 4*a*b.^2*(0+c2).*(x+c1).*(3+(0+c2).^2 + a.^2*(x+c1).^2) + b*(-3+(0+c2).^4 + 6*a.^2*(x+c1).^2 + a.^4*(x+c1).^4 + 2*(0+c2).^2.*(3 + a.^2*(x+c1).^2)) + 2*b.^3*(-1 + (0+c2).^4 + a.^2*(x+c1).^2 + (0+c2).^2*(4 + 3*a.^2*(x+c1).^2)));
dq_init = 1 + (0 + c2).^2 + b.^2 + (b*(0 + c2) + a*(x + c1)).^2;
q_init = (nq_init ./ (dq_init .^ 2));
ns_init = -(-3 + (0 + c2).^2 + b.^2*(-3 + (0 + c2).^2) + 4*1i*a*(x+c1) – 2*a*b*(0 + c2).*(x+c1) + a.^2*(x+c1).^2).*(exp(1i*(a*x + b*0)));
s_init = ns_init ./ dq_init;

initial_conditions = [q_init(:); s_init(:)];

t_span = T;

% Solve the differential equation numerically
options = odeset(‘RelTol’, 1e-4, ‘AbsTol’, 1e-8);
[t, u] = ode45(@(t, u) dydt(t, u, N), t_span, initial_conditions, options);

% Extract q and s from the solution
q = reshape(u(:, 1:N), [length(t), N]);
s = reshape(u(:, N+1:end), [length(t), N]);

% Calculate the exact solution over the mesh grid
[x, t] = meshgrid(X, T);
nq = -a*(b.^5*(1+(t + c2).^2).^2 – 8*a*(t+c2).*(x+c1) – 4*a*b.^4*(t+c2).*(1+(t+c2).^2).*(x+c1) – 4*a*b.^2*(t+c2).*(x+c1).*(3+(t+c2).^2 + a.^2*(x+c1).^2) + b*(-3+(t+c2).^4 + 6*a.^2*(x+c1).^2 + a.^4*(x+c1).^4 + 2*(t+c2).^2.*(3 + a.^2*(x+c1).^2)) + 2*b.^3*(-1 + (t+c2).^4 + a.^2*(x+c1).^2 + (t+c2).^2*(4 + 3*a.^2*(x+c1).^2)));
dq = 1 + (t + c2).^2 + b.^2 + (b*(t + c2) + a*(x + c1)).^2;
ns = -(-3 + (t + c2).^2 + b.^2*(-3 + (t + c2).^2) + 4*1i*a*(x+c1) – 2*a*b*(t + c2).*(x+c1) + a.^2*(x+c1).^2).*(exp(1i*(a*x + b*0)));
exact_q = (nq ./ (dq .^ 2));
exact_s = ns ./ dq;

% Compute the error
error_q = abs(q’ – exact_q);
error_s = abs(s’ – exact_s);

% Plotting the numerical solution
figure;
subplot(2, 1, 1);
surf(t, x, abs(q’));
title(‘Numerical Solution for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|q_{n+1} – q_n|’);
shading interp;

subplot(2, 1, 2);
surf(t, x, abs(s’));
title(‘Numerical Solution for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|s_{n+1} – s_n|’);
shading interp;

% 3D surface plot for exact solution
figure;
subplot(2, 1, 1);
surf(t, x, abs(exact_q’));
title(‘Exact Solution for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|q_{n+1} – q_n|’);
shading interp;

subplot(2, 1, 2);
surf(t, x, abs(exact_s’));
title(‘Exact Solution for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|s_{n+1} – s_n|’);
shading interp;

% Plot the error
figure;
subplot(2, 1, 1);
surf(t, x, error_q’);
title(‘Error for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘Error’);
shading interp;

subplot(2, 1, 2);
surf(t, x, error_s’);
title(‘Error for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘Error’);
shading interp;

% 2D error plot at specific times
figure;
times_to_plot = [0, 5, 10]; % Example times to plot error
for i = 1:length(times_to_plot)
time_idx = find(t_span >= times_to_plot(i), 1);
subplot(length(times_to_plot), 1, i);
plot(X, error_q(:, time_idx), ‘r’, X, error_s(:, time_idx), ‘b’);
title([‘Error at t = ‘, num2str(times_to_plot(i))]);
xlabel(‘x’);
ylabel(‘Error’);
legend(‘|q_{n+1} – q_n| Error’, ‘|s_{n+1} – s_n| Error’);
end

% Generate error table for fixed times from 0 to 10 in increments of 0.5
fixed_times = 0:0.5:10;
error_table_q = zeros(length(fixed_times), 1);
error_table_s = zeros(length(fixed_times), 1);

for i = 1:length(fixed_times)
time_idx = find(t_span >= fixed_times(i), 1);
error_table_q(i) = mean(error_q(:, time_idx));
error_table_s(i) = mean(error_s(:, time_idx));
end

% Display the error table
error_table = table(fixed_times’, error_table_q, error_table_s, …
‘VariableNames’, {‘Time’, ‘Error_q’, ‘Error_s’});
disp(error_table);

toc;

function dydt = dydt(~, u, N)
dydt = zeros(2*N, 1);
q = u(1:N);
s = u(N+1:2*N);

% Handle boundary conditions
for n = 3:N-2
dydt(n) = 1/2 * s(n) * conj(s(n)) – 1/2 * s(n+1) * conj(s(n+1));
dydt(N+n) = 1/2 * (q(n+1) – q(n)) * (s(n+1) – s(n));
end
endWhy numerical solutions is not consistent with the exact solutions? Also the error is too large. How to modify the code to get consistent results? See the attached image for the equation of motions and their exact solutions.
clear; clc;
tic;
c1 = 0; % real constant
c2 = 0; % real constant
a = 1; % real constant
b=0.1; % two values b=0.1, b=0.01
lambda = -1./(2*(1 + 1i)); % complex or imaginary
N = 1000;
X = linspace(-10, 10, N);
T = linspace(0, 10, N);

% Calculate initial conditions
[x, ~] = meshgrid(X, 0); % t = 0 for initial conditions
%psi_init = A * exp(-1i * lambda * x);
%phi_init = B * exp(1i * lambda * x);
nq_init = -a*(b.^5*(1+(0 + c2).^2).^2 – 8*a*(0+c2).*(x+c1) – 4*a*b.^4*(0+c2).*(1+(0+c2).^2).*(x+c1) – 4*a*b.^2*(0+c2).*(x+c1).*(3+(0+c2).^2 + a.^2*(x+c1).^2) + b*(-3+(0+c2).^4 + 6*a.^2*(x+c1).^2 + a.^4*(x+c1).^4 + 2*(0+c2).^2.*(3 + a.^2*(x+c1).^2)) + 2*b.^3*(-1 + (0+c2).^4 + a.^2*(x+c1).^2 + (0+c2).^2*(4 + 3*a.^2*(x+c1).^2)));
dq_init = 1 + (0 + c2).^2 + b.^2 + (b*(0 + c2) + a*(x + c1)).^2;
q_init = (nq_init ./ (dq_init .^ 2));
ns_init = -(-3 + (0 + c2).^2 + b.^2*(-3 + (0 + c2).^2) + 4*1i*a*(x+c1) – 2*a*b*(0 + c2).*(x+c1) + a.^2*(x+c1).^2).*(exp(1i*(a*x + b*0)));
s_init = ns_init ./ dq_init;

initial_conditions = [q_init(:); s_init(:)];

t_span = T;

% Solve the differential equation numerically
options = odeset(‘RelTol’, 1e-4, ‘AbsTol’, 1e-8);
[t, u] = ode45(@(t, u) dydt(t, u, N), t_span, initial_conditions, options);

% Extract q and s from the solution
q = reshape(u(:, 1:N), [length(t), N]);
s = reshape(u(:, N+1:end), [length(t), N]);

% Calculate the exact solution over the mesh grid
[x, t] = meshgrid(X, T);
nq = -a*(b.^5*(1+(t + c2).^2).^2 – 8*a*(t+c2).*(x+c1) – 4*a*b.^4*(t+c2).*(1+(t+c2).^2).*(x+c1) – 4*a*b.^2*(t+c2).*(x+c1).*(3+(t+c2).^2 + a.^2*(x+c1).^2) + b*(-3+(t+c2).^4 + 6*a.^2*(x+c1).^2 + a.^4*(x+c1).^4 + 2*(t+c2).^2.*(3 + a.^2*(x+c1).^2)) + 2*b.^3*(-1 + (t+c2).^4 + a.^2*(x+c1).^2 + (t+c2).^2*(4 + 3*a.^2*(x+c1).^2)));
dq = 1 + (t + c2).^2 + b.^2 + (b*(t + c2) + a*(x + c1)).^2;
ns = -(-3 + (t + c2).^2 + b.^2*(-3 + (t + c2).^2) + 4*1i*a*(x+c1) – 2*a*b*(t + c2).*(x+c1) + a.^2*(x+c1).^2).*(exp(1i*(a*x + b*0)));
exact_q = (nq ./ (dq .^ 2));
exact_s = ns ./ dq;

% Compute the error
error_q = abs(q’ – exact_q);
error_s = abs(s’ – exact_s);

% Plotting the numerical solution
figure;
subplot(2, 1, 1);
surf(t, x, abs(q’));
title(‘Numerical Solution for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|q_{n+1} – q_n|’);
shading interp;

subplot(2, 1, 2);
surf(t, x, abs(s’));
title(‘Numerical Solution for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|s_{n+1} – s_n|’);
shading interp;

% 3D surface plot for exact solution
figure;
subplot(2, 1, 1);
surf(t, x, abs(exact_q’));
title(‘Exact Solution for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|q_{n+1} – q_n|’);
shading interp;

subplot(2, 1, 2);
surf(t, x, abs(exact_s’));
title(‘Exact Solution for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|s_{n+1} – s_n|’);
shading interp;

% Plot the error
figure;
subplot(2, 1, 1);
surf(t, x, error_q’);
title(‘Error for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘Error’);
shading interp;

subplot(2, 1, 2);
surf(t, x, error_s’);
title(‘Error for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘Error’);
shading interp;

% 2D error plot at specific times
figure;
times_to_plot = [0, 5, 10]; % Example times to plot error
for i = 1:length(times_to_plot)
time_idx = find(t_span >= times_to_plot(i), 1);
subplot(length(times_to_plot), 1, i);
plot(X, error_q(:, time_idx), ‘r’, X, error_s(:, time_idx), ‘b’);
title([‘Error at t = ‘, num2str(times_to_plot(i))]);
xlabel(‘x’);
ylabel(‘Error’);
legend(‘|q_{n+1} – q_n| Error’, ‘|s_{n+1} – s_n| Error’);
end

% Generate error table for fixed times from 0 to 10 in increments of 0.5
fixed_times = 0:0.5:10;
error_table_q = zeros(length(fixed_times), 1);
error_table_s = zeros(length(fixed_times), 1);

for i = 1:length(fixed_times)
time_idx = find(t_span >= fixed_times(i), 1);
error_table_q(i) = mean(error_q(:, time_idx));
error_table_s(i) = mean(error_s(:, time_idx));
end

% Display the error table
error_table = table(fixed_times’, error_table_q, error_table_s, …
‘VariableNames’, {‘Time’, ‘Error_q’, ‘Error_s’});
disp(error_table);

toc;

function dydt = dydt(~, u, N)
dydt = zeros(2*N, 1);
q = u(1:N);
s = u(N+1:2*N);

% Handle boundary conditions
for n = 3:N-2
dydt(n) = 1/2 * s(n) * conj(s(n)) – 1/2 * s(n+1) * conj(s(n+1));
dydt(N+n) = 1/2 * (q(n+1) – q(n)) * (s(n+1) – s(n));
end
end Why numerical solutions is not consistent with the exact solutions? Also the error is too large. How to modify the code to get consistent results? See the attached image for the equation of motions and their exact solutions.
clear; clc;
tic;
c1 = 0; % real constant
c2 = 0; % real constant
a = 1; % real constant
b=0.1; % two values b=0.1, b=0.01
lambda = -1./(2*(1 + 1i)); % complex or imaginary
N = 1000;
X = linspace(-10, 10, N);
T = linspace(0, 10, N);

% Calculate initial conditions
[x, ~] = meshgrid(X, 0); % t = 0 for initial conditions
%psi_init = A * exp(-1i * lambda * x);
%phi_init = B * exp(1i * lambda * x);
nq_init = -a*(b.^5*(1+(0 + c2).^2).^2 – 8*a*(0+c2).*(x+c1) – 4*a*b.^4*(0+c2).*(1+(0+c2).^2).*(x+c1) – 4*a*b.^2*(0+c2).*(x+c1).*(3+(0+c2).^2 + a.^2*(x+c1).^2) + b*(-3+(0+c2).^4 + 6*a.^2*(x+c1).^2 + a.^4*(x+c1).^4 + 2*(0+c2).^2.*(3 + a.^2*(x+c1).^2)) + 2*b.^3*(-1 + (0+c2).^4 + a.^2*(x+c1).^2 + (0+c2).^2*(4 + 3*a.^2*(x+c1).^2)));
dq_init = 1 + (0 + c2).^2 + b.^2 + (b*(0 + c2) + a*(x + c1)).^2;
q_init = (nq_init ./ (dq_init .^ 2));
ns_init = -(-3 + (0 + c2).^2 + b.^2*(-3 + (0 + c2).^2) + 4*1i*a*(x+c1) – 2*a*b*(0 + c2).*(x+c1) + a.^2*(x+c1).^2).*(exp(1i*(a*x + b*0)));
s_init = ns_init ./ dq_init;

initial_conditions = [q_init(:); s_init(:)];

t_span = T;

% Solve the differential equation numerically
options = odeset(‘RelTol’, 1e-4, ‘AbsTol’, 1e-8);
[t, u] = ode45(@(t, u) dydt(t, u, N), t_span, initial_conditions, options);

% Extract q and s from the solution
q = reshape(u(:, 1:N), [length(t), N]);
s = reshape(u(:, N+1:end), [length(t), N]);

% Calculate the exact solution over the mesh grid
[x, t] = meshgrid(X, T);
nq = -a*(b.^5*(1+(t + c2).^2).^2 – 8*a*(t+c2).*(x+c1) – 4*a*b.^4*(t+c2).*(1+(t+c2).^2).*(x+c1) – 4*a*b.^2*(t+c2).*(x+c1).*(3+(t+c2).^2 + a.^2*(x+c1).^2) + b*(-3+(t+c2).^4 + 6*a.^2*(x+c1).^2 + a.^4*(x+c1).^4 + 2*(t+c2).^2.*(3 + a.^2*(x+c1).^2)) + 2*b.^3*(-1 + (t+c2).^4 + a.^2*(x+c1).^2 + (t+c2).^2*(4 + 3*a.^2*(x+c1).^2)));
dq = 1 + (t + c2).^2 + b.^2 + (b*(t + c2) + a*(x + c1)).^2;
ns = -(-3 + (t + c2).^2 + b.^2*(-3 + (t + c2).^2) + 4*1i*a*(x+c1) – 2*a*b*(t + c2).*(x+c1) + a.^2*(x+c1).^2).*(exp(1i*(a*x + b*0)));
exact_q = (nq ./ (dq .^ 2));
exact_s = ns ./ dq;

% Compute the error
error_q = abs(q’ – exact_q);
error_s = abs(s’ – exact_s);

% Plotting the numerical solution
figure;
subplot(2, 1, 1);
surf(t, x, abs(q’));
title(‘Numerical Solution for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|q_{n+1} – q_n|’);
shading interp;

subplot(2, 1, 2);
surf(t, x, abs(s’));
title(‘Numerical Solution for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|s_{n+1} – s_n|’);
shading interp;

% 3D surface plot for exact solution
figure;
subplot(2, 1, 1);
surf(t, x, abs(exact_q’));
title(‘Exact Solution for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|q_{n+1} – q_n|’);
shading interp;

subplot(2, 1, 2);
surf(t, x, abs(exact_s’));
title(‘Exact Solution for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘|s_{n+1} – s_n|’);
shading interp;

% Plot the error
figure;
subplot(2, 1, 1);
surf(t, x, error_q’);
title(‘Error for |q_{n+1} – q_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘Error’);
shading interp;

subplot(2, 1, 2);
surf(t, x, error_s’);
title(‘Error for |s_{n+1} – s_n|’);
xlabel(‘t’);
ylabel(‘x’);
zlabel(‘Error’);
shading interp;

% 2D error plot at specific times
figure;
times_to_plot = [0, 5, 10]; % Example times to plot error
for i = 1:length(times_to_plot)
time_idx = find(t_span >= times_to_plot(i), 1);
subplot(length(times_to_plot), 1, i);
plot(X, error_q(:, time_idx), ‘r’, X, error_s(:, time_idx), ‘b’);
title([‘Error at t = ‘, num2str(times_to_plot(i))]);
xlabel(‘x’);
ylabel(‘Error’);
legend(‘|q_{n+1} – q_n| Error’, ‘|s_{n+1} – s_n| Error’);
end

% Generate error table for fixed times from 0 to 10 in increments of 0.5
fixed_times = 0:0.5:10;
error_table_q = zeros(length(fixed_times), 1);
error_table_s = zeros(length(fixed_times), 1);

for i = 1:length(fixed_times)
time_idx = find(t_span >= fixed_times(i), 1);
error_table_q(i) = mean(error_q(:, time_idx));
error_table_s(i) = mean(error_s(:, time_idx));
end

% Display the error table
error_table = table(fixed_times’, error_table_q, error_table_s, …
‘VariableNames’, {‘Time’, ‘Error_q’, ‘Error_s’});
disp(error_table);

toc;

function dydt = dydt(~, u, N)
dydt = zeros(2*N, 1);
q = u(1:N);
s = u(N+1:2*N);

% Handle boundary conditions
for n = 3:N-2
dydt(n) = 1/2 * s(n) * conj(s(n)) – 1/2 * s(n+1) * conj(s(n+1));
dydt(N+n) = 1/2 * (q(n+1) – q(n)) * (s(n+1) – s(n));
end
end ode45, numerical integration, differential equations MATLAB Answers — New Questions

​