Category: Matlab
Category Archives: Matlab
The indices for validation and test sets are not being assigned correctly by ‘divideind’
I am new in matlab and I am trying to use ‘divideind’ to divide my data in order, but unfortunately it gives me this:
trainInd: [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 … ] (1×1418 double)
valInd: [1×0 double]
testInd: [1×0 double]
I would appreciate it if you tell me where is the problem?
this is my code:
X = data.inputs;
% Define the dependent variable
y = data.output;
Q = length(y);
% Define the split points
train_end = 1418;
val_end = 1773;
test_end = Q;
% Assign the indices
train_indices = 1:train_end;
val_indices = (train_end + 1):val_end;
test_indices = (val_end + 1):test_end;
% Use divideind to ensure proper division
[train_indices,val_indices,test_indices] = divideind(Q,train_indices,val_indices,test_indices);
% Divide the data into training and testing sets
x_train = log(X(train_indices, :));
x_valid = log(X(val_indices, :));
x_test = log(X(test_indices, :));
y_train = log(y(train_indices));
y_valid = log(y(val_indices));
y_test = log(y(test_indices));
% Transpose the data to fit the input format of the neural network
inputs_train = x_train’;
outputs_train = y_train’;
inputs_validation = x_valid’ ;
outputs_valid = y_valid’;
inputs_test = x_test’;
outputs_test = y_test’;
% Initialize arrays to store RMSE values
train_rmse = zeros(1, 10);
Valid_rmse = zeros(1, 10);
test_rmse = zeros(1, 10);
for i = 1:10
hiddenLayerSize = [i, i];
net = fitnet(hiddenLayerSize, ‘trainlm’);
net.trainParam.epochs = 200;
net.layers{1}.transferFcn = ‘poslin’;
net.layers{2}.transferFcn = ‘logsig’;
net.divideFcn = ‘divideind’;
net.divideMode = ‘sample’;
net.divideParam.trainInd = train_indices;
net.divideParam.valInd = val_indices;
net.divideParam.testInd = test_indices;
[net,tr] = train(net, inputs_train, outputs_train);
% Predict on training data
train_predictions = exp(net(inputs_train(:,tr.trainInd)));
yTrainTrue = exp(outputs_train(:,tr.trainInd));
train_rmse(i) = sqrt(mean((train_predictions – yTrainTrue).^2));
% Predict on Validation data
Valid_predictions = exp(net(inputs_validation(:,tr.valInd)));
yValTrue = exp(outputs_valid(:,tr.valInd));
Valid_rmse(i) = sqrt(mean((Valid_predictions – yValTrue).^2));
% Predict on Validation data
test_predictions = exp(net(inputs_test(:,tr.testInd)));
yTestTrue = exp(outputs_test(:,tr.testInd));
test_rmse(i) = sqrt(mean((test_predictions – yTestTrue).^2));
endI am new in matlab and I am trying to use ‘divideind’ to divide my data in order, but unfortunately it gives me this:
trainInd: [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 … ] (1×1418 double)
valInd: [1×0 double]
testInd: [1×0 double]
I would appreciate it if you tell me where is the problem?
this is my code:
X = data.inputs;
% Define the dependent variable
y = data.output;
Q = length(y);
% Define the split points
train_end = 1418;
val_end = 1773;
test_end = Q;
% Assign the indices
train_indices = 1:train_end;
val_indices = (train_end + 1):val_end;
test_indices = (val_end + 1):test_end;
% Use divideind to ensure proper division
[train_indices,val_indices,test_indices] = divideind(Q,train_indices,val_indices,test_indices);
% Divide the data into training and testing sets
x_train = log(X(train_indices, :));
x_valid = log(X(val_indices, :));
x_test = log(X(test_indices, :));
y_train = log(y(train_indices));
y_valid = log(y(val_indices));
y_test = log(y(test_indices));
% Transpose the data to fit the input format of the neural network
inputs_train = x_train’;
outputs_train = y_train’;
inputs_validation = x_valid’ ;
outputs_valid = y_valid’;
inputs_test = x_test’;
outputs_test = y_test’;
% Initialize arrays to store RMSE values
train_rmse = zeros(1, 10);
Valid_rmse = zeros(1, 10);
test_rmse = zeros(1, 10);
for i = 1:10
hiddenLayerSize = [i, i];
net = fitnet(hiddenLayerSize, ‘trainlm’);
net.trainParam.epochs = 200;
net.layers{1}.transferFcn = ‘poslin’;
net.layers{2}.transferFcn = ‘logsig’;
net.divideFcn = ‘divideind’;
net.divideMode = ‘sample’;
net.divideParam.trainInd = train_indices;
net.divideParam.valInd = val_indices;
net.divideParam.testInd = test_indices;
[net,tr] = train(net, inputs_train, outputs_train);
% Predict on training data
train_predictions = exp(net(inputs_train(:,tr.trainInd)));
yTrainTrue = exp(outputs_train(:,tr.trainInd));
train_rmse(i) = sqrt(mean((train_predictions – yTrainTrue).^2));
% Predict on Validation data
Valid_predictions = exp(net(inputs_validation(:,tr.valInd)));
yValTrue = exp(outputs_valid(:,tr.valInd));
Valid_rmse(i) = sqrt(mean((Valid_predictions – yValTrue).^2));
% Predict on Validation data
test_predictions = exp(net(inputs_test(:,tr.testInd)));
yTestTrue = exp(outputs_test(:,tr.testInd));
test_rmse(i) = sqrt(mean((test_predictions – yTestTrue).^2));
end I am new in matlab and I am trying to use ‘divideind’ to divide my data in order, but unfortunately it gives me this:
trainInd: [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 … ] (1×1418 double)
valInd: [1×0 double]
testInd: [1×0 double]
I would appreciate it if you tell me where is the problem?
this is my code:
X = data.inputs;
% Define the dependent variable
y = data.output;
Q = length(y);
% Define the split points
train_end = 1418;
val_end = 1773;
test_end = Q;
% Assign the indices
train_indices = 1:train_end;
val_indices = (train_end + 1):val_end;
test_indices = (val_end + 1):test_end;
% Use divideind to ensure proper division
[train_indices,val_indices,test_indices] = divideind(Q,train_indices,val_indices,test_indices);
% Divide the data into training and testing sets
x_train = log(X(train_indices, :));
x_valid = log(X(val_indices, :));
x_test = log(X(test_indices, :));
y_train = log(y(train_indices));
y_valid = log(y(val_indices));
y_test = log(y(test_indices));
% Transpose the data to fit the input format of the neural network
inputs_train = x_train’;
outputs_train = y_train’;
inputs_validation = x_valid’ ;
outputs_valid = y_valid’;
inputs_test = x_test’;
outputs_test = y_test’;
% Initialize arrays to store RMSE values
train_rmse = zeros(1, 10);
Valid_rmse = zeros(1, 10);
test_rmse = zeros(1, 10);
for i = 1:10
hiddenLayerSize = [i, i];
net = fitnet(hiddenLayerSize, ‘trainlm’);
net.trainParam.epochs = 200;
net.layers{1}.transferFcn = ‘poslin’;
net.layers{2}.transferFcn = ‘logsig’;
net.divideFcn = ‘divideind’;
net.divideMode = ‘sample’;
net.divideParam.trainInd = train_indices;
net.divideParam.valInd = val_indices;
net.divideParam.testInd = test_indices;
[net,tr] = train(net, inputs_train, outputs_train);
% Predict on training data
train_predictions = exp(net(inputs_train(:,tr.trainInd)));
yTrainTrue = exp(outputs_train(:,tr.trainInd));
train_rmse(i) = sqrt(mean((train_predictions – yTrainTrue).^2));
% Predict on Validation data
Valid_predictions = exp(net(inputs_validation(:,tr.valInd)));
yValTrue = exp(outputs_valid(:,tr.valInd));
Valid_rmse(i) = sqrt(mean((Valid_predictions – yValTrue).^2));
% Predict on Validation data
test_predictions = exp(net(inputs_test(:,tr.testInd)));
yTestTrue = exp(outputs_test(:,tr.testInd));
test_rmse(i) = sqrt(mean((test_predictions – yTestTrue).^2));
end divideind, ann, split data MATLAB Answers — New Questions
how to get MTF from PSF
Hi guys,
Can you please tell me how can I get MTF from the PSF. PSF is the matrix of 400X400 taken from the camera. Please see attached mat file.
Thank you very much.Hi guys,
Can you please tell me how can I get MTF from the PSF. PSF is the matrix of 400X400 taken from the camera. Please see attached mat file.
Thank you very much. Hi guys,
Can you please tell me how can I get MTF from the PSF. PSF is the matrix of 400X400 taken from the camera. Please see attached mat file.
Thank you very much. psf, mtf MATLAB Answers — New Questions
Are there any examples of using Transformers for Time Series?
Good afternoon All,
I’m researching how to develop/train/use Transformers (those similar to GPT) for S&P 500 historical data as an example. It seems like Transformers are mainly used for NLP but there are articles out there that utilize them for time series. Can any fellow lead the way in how could I get started in doing these in MATLAB? I have done general things like convolutions, LSTMs, etc. but not on a transformers caliber with the DL Toolbox.
I found this MATLAB package (pasted below) to be somewhat useful but it seems to lean more into language rather than historical data points.
https://www.mathworks.com/matlabcentral/fileexchange/107375-transformer-models
Another side question: It’s possibe to use the Deep Network Designer app to layout a transformers architecture?
Any tip is appreciated.
Thank you,
YerayGood afternoon All,
I’m researching how to develop/train/use Transformers (those similar to GPT) for S&P 500 historical data as an example. It seems like Transformers are mainly used for NLP but there are articles out there that utilize them for time series. Can any fellow lead the way in how could I get started in doing these in MATLAB? I have done general things like convolutions, LSTMs, etc. but not on a transformers caliber with the DL Toolbox.
I found this MATLAB package (pasted below) to be somewhat useful but it seems to lean more into language rather than historical data points.
https://www.mathworks.com/matlabcentral/fileexchange/107375-transformer-models
Another side question: It’s possibe to use the Deep Network Designer app to layout a transformers architecture?
Any tip is appreciated.
Thank you,
Yeray Good afternoon All,
I’m researching how to develop/train/use Transformers (those similar to GPT) for S&P 500 historical data as an example. It seems like Transformers are mainly used for NLP but there are articles out there that utilize them for time series. Can any fellow lead the way in how could I get started in doing these in MATLAB? I have done general things like convolutions, LSTMs, etc. but not on a transformers caliber with the DL Toolbox.
I found this MATLAB package (pasted below) to be somewhat useful but it seems to lean more into language rather than historical data points.
https://www.mathworks.com/matlabcentral/fileexchange/107375-transformer-models
Another side question: It’s possibe to use the Deep Network Designer app to layout a transformers architecture?
Any tip is appreciated.
Thank you,
Yeray deep learning, ai, transformers, gpt MATLAB Answers — New Questions
Problem with polyfit command (R2015a)
Below is a small script for using the polyfit command but surprisingly the last command gives me a completely wrong polynomial p. I don’t understand why. Thanks in advance for your help.
Below the script is the response from my version of MATLAB (R2015a).
%————————————————————————————————————–
%butta_sto_test
%
clear all
clc
x=[1:9]
y=[5,6,10,20,28,33,34,36,42]
p=polyfit(x,y,1)
[p,S]=polyfit(x,y,1)
[p,S,mu]=polyfit(x,y,1)
%———————————————————————————————————————-
x =
1 2 3 4 5 6 7 8 9
y =
5 6 10 20 28 33 34 36 42
p =
4.9833 -1.1389
p =
4.9833 -1.1389
S =
R: [2×2 double]
df: 7
normr: 8.4581
p =
13.6474 23.7778
S =
R: [2×2 double]
df: 7
normr: 8.4581
mu =
5.0000
2.7386Below is a small script for using the polyfit command but surprisingly the last command gives me a completely wrong polynomial p. I don’t understand why. Thanks in advance for your help.
Below the script is the response from my version of MATLAB (R2015a).
%————————————————————————————————————–
%butta_sto_test
%
clear all
clc
x=[1:9]
y=[5,6,10,20,28,33,34,36,42]
p=polyfit(x,y,1)
[p,S]=polyfit(x,y,1)
[p,S,mu]=polyfit(x,y,1)
%———————————————————————————————————————-
x =
1 2 3 4 5 6 7 8 9
y =
5 6 10 20 28 33 34 36 42
p =
4.9833 -1.1389
p =
4.9833 -1.1389
S =
R: [2×2 double]
df: 7
normr: 8.4581
p =
13.6474 23.7778
S =
R: [2×2 double]
df: 7
normr: 8.4581
mu =
5.0000
2.7386 Below is a small script for using the polyfit command but surprisingly the last command gives me a completely wrong polynomial p. I don’t understand why. Thanks in advance for your help.
Below the script is the response from my version of MATLAB (R2015a).
%————————————————————————————————————–
%butta_sto_test
%
clear all
clc
x=[1:9]
y=[5,6,10,20,28,33,34,36,42]
p=polyfit(x,y,1)
[p,S]=polyfit(x,y,1)
[p,S,mu]=polyfit(x,y,1)
%———————————————————————————————————————-
x =
1 2 3 4 5 6 7 8 9
y =
5 6 10 20 28 33 34 36 42
p =
4.9833 -1.1389
p =
4.9833 -1.1389
S =
R: [2×2 double]
df: 7
normr: 8.4581
p =
13.6474 23.7778
S =
R: [2×2 double]
df: 7
normr: 8.4581
mu =
5.0000
2.7386 polyfit MATLAB Answers — New Questions
I have a PSF code and I need to generate the MTF and plot it along with the spatial frequency graph, similar to the image attached. I would greatly appreciate your assistance
% Imaging with Zone Plate
for tt=1:1
z2 = s2+s2/5000*tt;
h2 = exp(1i*2*pi*z2*sqrt(1/(lambda)^2-(x/L).^2-(y/L).^2));
E2 = ifft2(fftshift(fftshift(fft2(E.*exp(1i*FZP))).*h2));
psf=fftshift(fft2(exp(1i*FZP)));
resulting_image = abs(E2).^2;
figure(200+tt)
imagesc(abs(E2).^2);
title(‘resulting_image’)
axis square;
I need assistance in calculating MTFfrom PSF, plotting MTF vs spatial frequency graph, and calculating the Strehl ratio. I am a beginner and have been struggling with this for a long time. Your help would be greatly appreciated.% Imaging with Zone Plate
for tt=1:1
z2 = s2+s2/5000*tt;
h2 = exp(1i*2*pi*z2*sqrt(1/(lambda)^2-(x/L).^2-(y/L).^2));
E2 = ifft2(fftshift(fftshift(fft2(E.*exp(1i*FZP))).*h2));
psf=fftshift(fft2(exp(1i*FZP)));
resulting_image = abs(E2).^2;
figure(200+tt)
imagesc(abs(E2).^2);
title(‘resulting_image’)
axis square;
I need assistance in calculating MTFfrom PSF, plotting MTF vs spatial frequency graph, and calculating the Strehl ratio. I am a beginner and have been struggling with this for a long time. Your help would be greatly appreciated. % Imaging with Zone Plate
for tt=1:1
z2 = s2+s2/5000*tt;
h2 = exp(1i*2*pi*z2*sqrt(1/(lambda)^2-(x/L).^2-(y/L).^2));
E2 = ifft2(fftshift(fftshift(fft2(E.*exp(1i*FZP))).*h2));
psf=fftshift(fft2(exp(1i*FZP)));
resulting_image = abs(E2).^2;
figure(200+tt)
imagesc(abs(E2).^2);
title(‘resulting_image’)
axis square;
I need assistance in calculating MTFfrom PSF, plotting MTF vs spatial frequency graph, and calculating the Strehl ratio. I am a beginner and have been struggling with this for a long time. Your help would be greatly appreciated. mtf, psf MATLAB Answers — New Questions
How can I customize a built-in function “HelperOrientationViewer()”
I’m trying to get the real-time orientation of the IMU9250 sensor by following the documentation: "Estimating Orientation Using Inertial Sensor Fusion and MPU-9250"
Can I customize the HelperOrientationViewer() function? I want to add another 3D object to the 3D plot. I also want to change the normal 3D object to another mesh STL object. and I want to change the scale of each axis to be more bigger. How can I do these?
I tryed to use poseplot() function, but it is too slow comparing to HelperOrientationViewer() functionI’m trying to get the real-time orientation of the IMU9250 sensor by following the documentation: "Estimating Orientation Using Inertial Sensor Fusion and MPU-9250"
Can I customize the HelperOrientationViewer() function? I want to add another 3D object to the 3D plot. I also want to change the normal 3D object to another mesh STL object. and I want to change the scale of each axis to be more bigger. How can I do these?
I tryed to use poseplot() function, but it is too slow comparing to HelperOrientationViewer() function I’m trying to get the real-time orientation of the IMU9250 sensor by following the documentation: "Estimating Orientation Using Inertial Sensor Fusion and MPU-9250"
Can I customize the HelperOrientationViewer() function? I want to add another 3D object to the 3D plot. I also want to change the normal 3D object to another mesh STL object. and I want to change the scale of each axis to be more bigger. How can I do these?
I tryed to use poseplot() function, but it is too slow comparing to HelperOrientationViewer() function arduino, sensor fusion, helperorientationviewer, imu9250, graphics, poseplot MATLAB Answers — New Questions
How to stop quadratic formula calculator from giving inverse outputs?
Hi, I’m doing a worksheet based around matlab which involved solving quadratic formulas. I probably overcomplicated it and now made a calculator that keeps giving the inverse of my expected answer, could anyone help please?
for k = 1 : 7
disp ( ‘ For the equation ax^2 + bx + c ‘ )
a = input ( ‘ Enter a; ‘ );
b = input ( ‘ Enter b; ‘ );
c = input ( ‘ Enter c; ‘ );
D = ( b^2 ) – ( 4 * a * c );
if D < 0
fprintf ( ‘ n The equation has no real roots . nn ‘)
elseif D == 0
root = -b / ( 2 * a );
fprintf ( ‘ nThe equation has one root, n ‘ )
fprintf ( ‘ %.7fnn ‘ , root )
else
r1 = ( – b + sqrt ( D ) ) / ( 2 * a ) ;
r2 = ( – b – sqrt ( D ) ) / ( 2 * a ) ;
fprintf ( ‘n The equation has two roots, n ‘ )
fprintf ( ‘ %.7f and %.7fnn ‘ , r1 ,r2 )
end
end
For the equation ax^2 + bx + c
Enter a;
1
Enter b;
9
Enter c;
14
The equation has two roots,
-2.0000000 and -7.0000000
For the equation ax^2 + bx + c
Enter a;
1
Enter b;
-3
Enter c;
-18
The equation has two roots,
6.0000000 and -3.0000000Hi, I’m doing a worksheet based around matlab which involved solving quadratic formulas. I probably overcomplicated it and now made a calculator that keeps giving the inverse of my expected answer, could anyone help please?
for k = 1 : 7
disp ( ‘ For the equation ax^2 + bx + c ‘ )
a = input ( ‘ Enter a; ‘ );
b = input ( ‘ Enter b; ‘ );
c = input ( ‘ Enter c; ‘ );
D = ( b^2 ) – ( 4 * a * c );
if D < 0
fprintf ( ‘ n The equation has no real roots . nn ‘)
elseif D == 0
root = -b / ( 2 * a );
fprintf ( ‘ nThe equation has one root, n ‘ )
fprintf ( ‘ %.7fnn ‘ , root )
else
r1 = ( – b + sqrt ( D ) ) / ( 2 * a ) ;
r2 = ( – b – sqrt ( D ) ) / ( 2 * a ) ;
fprintf ( ‘n The equation has two roots, n ‘ )
fprintf ( ‘ %.7f and %.7fnn ‘ , r1 ,r2 )
end
end
For the equation ax^2 + bx + c
Enter a;
1
Enter b;
9
Enter c;
14
The equation has two roots,
-2.0000000 and -7.0000000
For the equation ax^2 + bx + c
Enter a;
1
Enter b;
-3
Enter c;
-18
The equation has two roots,
6.0000000 and -3.0000000 Hi, I’m doing a worksheet based around matlab which involved solving quadratic formulas. I probably overcomplicated it and now made a calculator that keeps giving the inverse of my expected answer, could anyone help please?
for k = 1 : 7
disp ( ‘ For the equation ax^2 + bx + c ‘ )
a = input ( ‘ Enter a; ‘ );
b = input ( ‘ Enter b; ‘ );
c = input ( ‘ Enter c; ‘ );
D = ( b^2 ) – ( 4 * a * c );
if D < 0
fprintf ( ‘ n The equation has no real roots . nn ‘)
elseif D == 0
root = -b / ( 2 * a );
fprintf ( ‘ nThe equation has one root, n ‘ )
fprintf ( ‘ %.7fnn ‘ , root )
else
r1 = ( – b + sqrt ( D ) ) / ( 2 * a ) ;
r2 = ( – b – sqrt ( D ) ) / ( 2 * a ) ;
fprintf ( ‘n The equation has two roots, n ‘ )
fprintf ( ‘ %.7f and %.7fnn ‘ , r1 ,r2 )
end
end
For the equation ax^2 + bx + c
Enter a;
1
Enter b;
9
Enter c;
14
The equation has two roots,
-2.0000000 and -7.0000000
For the equation ax^2 + bx + c
Enter a;
1
Enter b;
-3
Enter c;
-18
The equation has two roots,
6.0000000 and -3.0000000 quadratic formula MATLAB Answers — New Questions
Write a function called freezing that takes a vector of numbers that correspond to daily low temperatures in Fahrenheit. Return numfreeze, the number of days with sub freezing
Please help!
3 is correct for Test 1, but I am going crazy trying to figure out why 5 is incorrect for Test 2.
function numfreeze = freezing(a)
numfreeze = sum(a>=32);
endPlease help!
3 is correct for Test 1, but I am going crazy trying to figure out why 5 is incorrect for Test 2.
function numfreeze = freezing(a)
numfreeze = sum(a>=32);
end Please help!
3 is correct for Test 1, but I am going crazy trying to figure out why 5 is incorrect for Test 2.
function numfreeze = freezing(a)
numfreeze = sum(a>=32);
end logical array MATLAB Answers — New Questions
Error when running function with batch()
Hi. I am using batch() to run a function with a par-for loop inside. When I run it as a script (instead of a function), everything works fine. However, when I replace the script with the equivalent function, I get an error. Here are the details:
%(A) Using batch() to call script. This works fine!
job = batch(parcluster(‘local’), ‘batchscript’,’Pool’,5);
%(B) Using batch() to call function. This produces an error.
jm=6;
job = batch(parcluster(‘local’),’batchfun’,0,{jm},’Pool’,5);
The script (batchscript.m) I am running in A, looks something like that:
jm = 6;
run(‘Paths.m’) %This creates ‘pth’ using ‘jm’
parfor vh = 1:length(alld)
parfun01(vh,jm,pth)
end
The function (batchfun.m) I am running in B, looks something like that:
function batchfun(jm)
run(‘Paths.m’) %This creates ‘pth’ using ‘jm’
parfor vh = 1:length(alld)
parfun01(vh,jm,pth)
end
They reason why I decided to turn this from a scrip to a function was because I wanted to pass jm as an input, instead of always changing the script. However, I cannot understand why it works perfectly fine as a script, but not as a function.
The error I am getting when running B) is the following:
getReport(job.Tasks(1).Error)
ans =
‘Error using batchfun (line 66)
The source code (E:Analysisbatchfun.m) for the parfor-loop that is trying to
execute on the worker could not be found.
Caused by:
Unrecognized function or variable ‘pth’.
Error using remoteParallelFunction (line 84)
Worker unable to find file.
Unrecognized function or variable ‘pth’.’
Anyone could please help me with what’s wrong? Thanks!Hi. I am using batch() to run a function with a par-for loop inside. When I run it as a script (instead of a function), everything works fine. However, when I replace the script with the equivalent function, I get an error. Here are the details:
%(A) Using batch() to call script. This works fine!
job = batch(parcluster(‘local’), ‘batchscript’,’Pool’,5);
%(B) Using batch() to call function. This produces an error.
jm=6;
job = batch(parcluster(‘local’),’batchfun’,0,{jm},’Pool’,5);
The script (batchscript.m) I am running in A, looks something like that:
jm = 6;
run(‘Paths.m’) %This creates ‘pth’ using ‘jm’
parfor vh = 1:length(alld)
parfun01(vh,jm,pth)
end
The function (batchfun.m) I am running in B, looks something like that:
function batchfun(jm)
run(‘Paths.m’) %This creates ‘pth’ using ‘jm’
parfor vh = 1:length(alld)
parfun01(vh,jm,pth)
end
They reason why I decided to turn this from a scrip to a function was because I wanted to pass jm as an input, instead of always changing the script. However, I cannot understand why it works perfectly fine as a script, but not as a function.
The error I am getting when running B) is the following:
getReport(job.Tasks(1).Error)
ans =
‘Error using batchfun (line 66)
The source code (E:Analysisbatchfun.m) for the parfor-loop that is trying to
execute on the worker could not be found.
Caused by:
Unrecognized function or variable ‘pth’.
Error using remoteParallelFunction (line 84)
Worker unable to find file.
Unrecognized function or variable ‘pth’.’
Anyone could please help me with what’s wrong? Thanks! Hi. I am using batch() to run a function with a par-for loop inside. When I run it as a script (instead of a function), everything works fine. However, when I replace the script with the equivalent function, I get an error. Here are the details:
%(A) Using batch() to call script. This works fine!
job = batch(parcluster(‘local’), ‘batchscript’,’Pool’,5);
%(B) Using batch() to call function. This produces an error.
jm=6;
job = batch(parcluster(‘local’),’batchfun’,0,{jm},’Pool’,5);
The script (batchscript.m) I am running in A, looks something like that:
jm = 6;
run(‘Paths.m’) %This creates ‘pth’ using ‘jm’
parfor vh = 1:length(alld)
parfun01(vh,jm,pth)
end
The function (batchfun.m) I am running in B, looks something like that:
function batchfun(jm)
run(‘Paths.m’) %This creates ‘pth’ using ‘jm’
parfor vh = 1:length(alld)
parfun01(vh,jm,pth)
end
They reason why I decided to turn this from a scrip to a function was because I wanted to pass jm as an input, instead of always changing the script. However, I cannot understand why it works perfectly fine as a script, but not as a function.
The error I am getting when running B) is the following:
getReport(job.Tasks(1).Error)
ans =
‘Error using batchfun (line 66)
The source code (E:Analysisbatchfun.m) for the parfor-loop that is trying to
execute on the worker could not be found.
Caused by:
Unrecognized function or variable ‘pth’.
Error using remoteParallelFunction (line 84)
Worker unable to find file.
Unrecognized function or variable ‘pth’.’
Anyone could please help me with what’s wrong? Thanks! batch, parfor, parallel computing, parallel MATLAB Answers — New Questions
using hooks with git source control from Matlab. Hooks that work from outside Matlab are not executed from within.
Hello Everyone,
I’m working with Matlab 2016a try to get the source control going. Commit, Push, Fetch etc work fine, from both the context menu in the Current Folder Window and from the GUI shipped with git.
But when I implement a server-sided hook in the remote repository, Matlab fails me by ignoring the hook. When I commit from outside Matlab, the hook works.
hook-file is located in ..remote_repohooks by the name ‘pre-receive’ without extension. Source Code:
#!/usr/bin/env python
import matplotlib.pyplot as plt
import sys
plt.figure()
plt.show()
sys.exit(1)
The hook is written for testing purposes only. It opens a matplotlib window and when the window is closed, an error is thrown as the Push-Operation is aborted.
As I mentioned, when I push my commit from outside Matlab via the GUI, the expected happens, matplotlib-window and error. But when I push from the same local repository to the same remote repository, the push is successful…
Is there a way to execute hooks both server-sided and local ones from Matlab?
Thanks in advance and greetings,
TobiasHello Everyone,
I’m working with Matlab 2016a try to get the source control going. Commit, Push, Fetch etc work fine, from both the context menu in the Current Folder Window and from the GUI shipped with git.
But when I implement a server-sided hook in the remote repository, Matlab fails me by ignoring the hook. When I commit from outside Matlab, the hook works.
hook-file is located in ..remote_repohooks by the name ‘pre-receive’ without extension. Source Code:
#!/usr/bin/env python
import matplotlib.pyplot as plt
import sys
plt.figure()
plt.show()
sys.exit(1)
The hook is written for testing purposes only. It opens a matplotlib window and when the window is closed, an error is thrown as the Push-Operation is aborted.
As I mentioned, when I push my commit from outside Matlab via the GUI, the expected happens, matplotlib-window and error. But when I push from the same local repository to the same remote repository, the push is successful…
Is there a way to execute hooks both server-sided and local ones from Matlab?
Thanks in advance and greetings,
Tobias Hello Everyone,
I’m working with Matlab 2016a try to get the source control going. Commit, Push, Fetch etc work fine, from both the context menu in the Current Folder Window and from the GUI shipped with git.
But when I implement a server-sided hook in the remote repository, Matlab fails me by ignoring the hook. When I commit from outside Matlab, the hook works.
hook-file is located in ..remote_repohooks by the name ‘pre-receive’ without extension. Source Code:
#!/usr/bin/env python
import matplotlib.pyplot as plt
import sys
plt.figure()
plt.show()
sys.exit(1)
The hook is written for testing purposes only. It opens a matplotlib window and when the window is closed, an error is thrown as the Push-Operation is aborted.
As I mentioned, when I push my commit from outside Matlab via the GUI, the expected happens, matplotlib-window and error. But when I push from the same local repository to the same remote repository, the push is successful…
Is there a way to execute hooks both server-sided and local ones from Matlab?
Thanks in advance and greetings,
Tobias source control, git, hooks MATLAB Answers — New Questions
data not generating and error in plotting
% Define initial parameters
lambda_init = 2;
kappa_init = 1;
theta_k_init = pi/10;
R_init = 7;
rout = 3;
% Define c
c = sqrt(4 – (rout/R_init)^2);
% Define the function for the differential equations
f = @(r, y) [y(2); ((-1*(lambda_init+1)*(r*y(2)+y(1)))+(1/r)*((kappa_init*r^2*cos(theta_k_init)-(lambda_init+1))*y(1))-(kappa_init*r*y(3)*sin(theta_k_init))+(-16*lambda_init*r^2)/(c^4*R_init^2)); y(4); ((-1*(r*y(4)+y(3)))+(1/r)*((kappa_init*r^2*cos(theta_k_init)-1)*y(3))+(kappa_init*r*y(1)*sin(theta_k_init)))];
% Solve the differential equations using ode45
r_span = linspace(0, rout, 100); % Define the range of r values
[~, sol] = ode45(f, r_span, [0, 0, 0, 0]);
% Extract solutions
dr_sol = sol(:,1);
dtheta_sol = sol(:,3);
% Plot the solutions
figure;
subplot(2,1,1);
plot(r_span, dr_sol, ‘b-‘);
xlabel(‘r’);
ylabel(‘dr(r)’);
title(‘Solution of dr(r) vs r’);
subplot(2,1,2);
plot(r_span, dtheta_sol, ‘b-‘);
xlabel(‘r’);
ylabel(‘dtheta(r)’);
title(‘Solution of dtheta(r) vs r’);
I am trying to solve the coupled differential eqns for d_r(r) vs r and d_theta(r) vs r for these parameter values and boundary conditions so that it is zero at both end d_r(0) = d_r(rout) = 0 and the same for d_theta. However, Not able to see the plot and data. please suggest me errors.% Define initial parameters
lambda_init = 2;
kappa_init = 1;
theta_k_init = pi/10;
R_init = 7;
rout = 3;
% Define c
c = sqrt(4 – (rout/R_init)^2);
% Define the function for the differential equations
f = @(r, y) [y(2); ((-1*(lambda_init+1)*(r*y(2)+y(1)))+(1/r)*((kappa_init*r^2*cos(theta_k_init)-(lambda_init+1))*y(1))-(kappa_init*r*y(3)*sin(theta_k_init))+(-16*lambda_init*r^2)/(c^4*R_init^2)); y(4); ((-1*(r*y(4)+y(3)))+(1/r)*((kappa_init*r^2*cos(theta_k_init)-1)*y(3))+(kappa_init*r*y(1)*sin(theta_k_init)))];
% Solve the differential equations using ode45
r_span = linspace(0, rout, 100); % Define the range of r values
[~, sol] = ode45(f, r_span, [0, 0, 0, 0]);
% Extract solutions
dr_sol = sol(:,1);
dtheta_sol = sol(:,3);
% Plot the solutions
figure;
subplot(2,1,1);
plot(r_span, dr_sol, ‘b-‘);
xlabel(‘r’);
ylabel(‘dr(r)’);
title(‘Solution of dr(r) vs r’);
subplot(2,1,2);
plot(r_span, dtheta_sol, ‘b-‘);
xlabel(‘r’);
ylabel(‘dtheta(r)’);
title(‘Solution of dtheta(r) vs r’);
I am trying to solve the coupled differential eqns for d_r(r) vs r and d_theta(r) vs r for these parameter values and boundary conditions so that it is zero at both end d_r(0) = d_r(rout) = 0 and the same for d_theta. However, Not able to see the plot and data. please suggest me errors. % Define initial parameters
lambda_init = 2;
kappa_init = 1;
theta_k_init = pi/10;
R_init = 7;
rout = 3;
% Define c
c = sqrt(4 – (rout/R_init)^2);
% Define the function for the differential equations
f = @(r, y) [y(2); ((-1*(lambda_init+1)*(r*y(2)+y(1)))+(1/r)*((kappa_init*r^2*cos(theta_k_init)-(lambda_init+1))*y(1))-(kappa_init*r*y(3)*sin(theta_k_init))+(-16*lambda_init*r^2)/(c^4*R_init^2)); y(4); ((-1*(r*y(4)+y(3)))+(1/r)*((kappa_init*r^2*cos(theta_k_init)-1)*y(3))+(kappa_init*r*y(1)*sin(theta_k_init)))];
% Solve the differential equations using ode45
r_span = linspace(0, rout, 100); % Define the range of r values
[~, sol] = ode45(f, r_span, [0, 0, 0, 0]);
% Extract solutions
dr_sol = sol(:,1);
dtheta_sol = sol(:,3);
% Plot the solutions
figure;
subplot(2,1,1);
plot(r_span, dr_sol, ‘b-‘);
xlabel(‘r’);
ylabel(‘dr(r)’);
title(‘Solution of dr(r) vs r’);
subplot(2,1,2);
plot(r_span, dtheta_sol, ‘b-‘);
xlabel(‘r’);
ylabel(‘dtheta(r)’);
title(‘Solution of dtheta(r) vs r’);
I am trying to solve the coupled differential eqns for d_r(r) vs r and d_theta(r) vs r for these parameter values and boundary conditions so that it is zero at both end d_r(0) = d_r(rout) = 0 and the same for d_theta. However, Not able to see the plot and data. please suggest me errors. eqn solve, plotting MATLAB Answers — New Questions
Manual and Automatic Image Recognition
I need to identify bright and dark spots in an image, but the contrast of these spots is not very good. Processing the image by taking the maximum and minimum contrast values from a portion of the image is always difficult. Someone told me that defining the circularity of contrast could help identify these spots, and this method is feasible. My code is provided below.
My question is: I want to manually designate a region, and then have the program automatically re-identify the bright and dark spots in the designated area after it finishes identifying the rest of the image. For example, in the image below, is there a way to help me identify the black and white circles? Additionally, I want to apply a Voronoi Diagram to the identified spots in the image.
code:
clc
clear
% read image
inpict = imread(‘image.png’);
inpict = im2gray(inpict);
% logfilter image
t=1;
sigma=sqrt(t);
Fsize = ceil(sigma*3)*2+1;
h1 = fspecial(‘gaussian’,[Fsize,Fsize],sigma);
h2 = fspecial(‘log’,[Fsize,Fsize],sigma);
filter1 = uint8(imfilter(inpict,h1,’replicate’));
filter2 = t * imfilter(inpict, h2, ‘replicate’);
bw = filter2 >=0;
filter2= uint8(filter2 + 128);
figure
subplot(221), imshow(filter1);
title(sprintf(‘filtered by Gaussinan, t =%d’,t));
subplot(222), imshow(filter2, []);
title(sprintf(‘filtered by log, t =%d’,t));
subplot(223), imshow(bw*255);
title(sprintf(‘binarized, t =%d’,t));
%Automatic Image Recognition
inpict = imflatfield(filter1,30); %Flat-field correction of images
inpict = adapthisteq(inpict);
inpict = imresize(inpict, 2);
inpict = imadjust(inpict);
%Automatic Image Recognition (1.Creating Local Extreme Masks)
minex = 30; %black
maxex = 20; %white
maskmx = imextendedmax(inpict, minex); % black
maskmn = imextendedmin(inpict, maxex); % white
% 2.Removal of part of the plaque on the border
maskmx = imclearborder(maskmx);
maskmn = imclearborder(maskmn);
% filtering based on patch attributes
minc = 0.2;
maxecc = 0.87;
minarea = 10;
% Maximum Mask
Lmx = bwlabel(maskmx);
Smx = regionprops(maskmx, ‘circularity’, ‘eccentricity’, ‘area’);
goodblobs = find([Smx.Circularity] >= minc …
& [Smx.Eccentricity] <= maxecc …
& [Smx.Area] >= minarea);
maskmx = any(Lmx == permute(goodblobs, [1 3 2]), 3);
% minimum value mask
Lmn = bwlabel(maskmn);
Smn = regionprops(maskmn, ‘circularity’, ‘eccentricity’, ‘area’);
goodblobs = find([Smn.Circularity] >= minc …
& [Smn.Eccentricity] <= maxecc …
& [Smn.Area] >= minarea);
maskmn = any(Lmn == permute(goodblobs, [1 3 2]), 3);
% Visualisation masks and images
outpict = cat(3, im2uint8(maskmx), im2uint8(maskmn), inpict);
imshow(outpict, ‘border’, ‘tight’);I need to identify bright and dark spots in an image, but the contrast of these spots is not very good. Processing the image by taking the maximum and minimum contrast values from a portion of the image is always difficult. Someone told me that defining the circularity of contrast could help identify these spots, and this method is feasible. My code is provided below.
My question is: I want to manually designate a region, and then have the program automatically re-identify the bright and dark spots in the designated area after it finishes identifying the rest of the image. For example, in the image below, is there a way to help me identify the black and white circles? Additionally, I want to apply a Voronoi Diagram to the identified spots in the image.
code:
clc
clear
% read image
inpict = imread(‘image.png’);
inpict = im2gray(inpict);
% logfilter image
t=1;
sigma=sqrt(t);
Fsize = ceil(sigma*3)*2+1;
h1 = fspecial(‘gaussian’,[Fsize,Fsize],sigma);
h2 = fspecial(‘log’,[Fsize,Fsize],sigma);
filter1 = uint8(imfilter(inpict,h1,’replicate’));
filter2 = t * imfilter(inpict, h2, ‘replicate’);
bw = filter2 >=0;
filter2= uint8(filter2 + 128);
figure
subplot(221), imshow(filter1);
title(sprintf(‘filtered by Gaussinan, t =%d’,t));
subplot(222), imshow(filter2, []);
title(sprintf(‘filtered by log, t =%d’,t));
subplot(223), imshow(bw*255);
title(sprintf(‘binarized, t =%d’,t));
%Automatic Image Recognition
inpict = imflatfield(filter1,30); %Flat-field correction of images
inpict = adapthisteq(inpict);
inpict = imresize(inpict, 2);
inpict = imadjust(inpict);
%Automatic Image Recognition (1.Creating Local Extreme Masks)
minex = 30; %black
maxex = 20; %white
maskmx = imextendedmax(inpict, minex); % black
maskmn = imextendedmin(inpict, maxex); % white
% 2.Removal of part of the plaque on the border
maskmx = imclearborder(maskmx);
maskmn = imclearborder(maskmn);
% filtering based on patch attributes
minc = 0.2;
maxecc = 0.87;
minarea = 10;
% Maximum Mask
Lmx = bwlabel(maskmx);
Smx = regionprops(maskmx, ‘circularity’, ‘eccentricity’, ‘area’);
goodblobs = find([Smx.Circularity] >= minc …
& [Smx.Eccentricity] <= maxecc …
& [Smx.Area] >= minarea);
maskmx = any(Lmx == permute(goodblobs, [1 3 2]), 3);
% minimum value mask
Lmn = bwlabel(maskmn);
Smn = regionprops(maskmn, ‘circularity’, ‘eccentricity’, ‘area’);
goodblobs = find([Smn.Circularity] >= minc …
& [Smn.Eccentricity] <= maxecc …
& [Smn.Area] >= minarea);
maskmn = any(Lmn == permute(goodblobs, [1 3 2]), 3);
% Visualisation masks and images
outpict = cat(3, im2uint8(maskmx), im2uint8(maskmn), inpict);
imshow(outpict, ‘border’, ‘tight’); I need to identify bright and dark spots in an image, but the contrast of these spots is not very good. Processing the image by taking the maximum and minimum contrast values from a portion of the image is always difficult. Someone told me that defining the circularity of contrast could help identify these spots, and this method is feasible. My code is provided below.
My question is: I want to manually designate a region, and then have the program automatically re-identify the bright and dark spots in the designated area after it finishes identifying the rest of the image. For example, in the image below, is there a way to help me identify the black and white circles? Additionally, I want to apply a Voronoi Diagram to the identified spots in the image.
code:
clc
clear
% read image
inpict = imread(‘image.png’);
inpict = im2gray(inpict);
% logfilter image
t=1;
sigma=sqrt(t);
Fsize = ceil(sigma*3)*2+1;
h1 = fspecial(‘gaussian’,[Fsize,Fsize],sigma);
h2 = fspecial(‘log’,[Fsize,Fsize],sigma);
filter1 = uint8(imfilter(inpict,h1,’replicate’));
filter2 = t * imfilter(inpict, h2, ‘replicate’);
bw = filter2 >=0;
filter2= uint8(filter2 + 128);
figure
subplot(221), imshow(filter1);
title(sprintf(‘filtered by Gaussinan, t =%d’,t));
subplot(222), imshow(filter2, []);
title(sprintf(‘filtered by log, t =%d’,t));
subplot(223), imshow(bw*255);
title(sprintf(‘binarized, t =%d’,t));
%Automatic Image Recognition
inpict = imflatfield(filter1,30); %Flat-field correction of images
inpict = adapthisteq(inpict);
inpict = imresize(inpict, 2);
inpict = imadjust(inpict);
%Automatic Image Recognition (1.Creating Local Extreme Masks)
minex = 30; %black
maxex = 20; %white
maskmx = imextendedmax(inpict, minex); % black
maskmn = imextendedmin(inpict, maxex); % white
% 2.Removal of part of the plaque on the border
maskmx = imclearborder(maskmx);
maskmn = imclearborder(maskmn);
% filtering based on patch attributes
minc = 0.2;
maxecc = 0.87;
minarea = 10;
% Maximum Mask
Lmx = bwlabel(maskmx);
Smx = regionprops(maskmx, ‘circularity’, ‘eccentricity’, ‘area’);
goodblobs = find([Smx.Circularity] >= minc …
& [Smx.Eccentricity] <= maxecc …
& [Smx.Area] >= minarea);
maskmx = any(Lmx == permute(goodblobs, [1 3 2]), 3);
% minimum value mask
Lmn = bwlabel(maskmn);
Smn = regionprops(maskmn, ‘circularity’, ‘eccentricity’, ‘area’);
goodblobs = find([Smn.Circularity] >= minc …
& [Smn.Eccentricity] <= maxecc …
& [Smn.Area] >= minarea);
maskmn = any(Lmn == permute(goodblobs, [1 3 2]), 3);
% Visualisation masks and images
outpict = cat(3, im2uint8(maskmx), im2uint8(maskmn), inpict);
imshow(outpict, ‘border’, ‘tight’); image analysis, recognition, graphics, image segmentation, skiz MATLAB Answers — New Questions
Second order coupled differential equation
I need help in solving the second order coupled differential equations as I’m unable to get the right solution after spending many days.
I want to solve it with two for loops (One for frequency and other for Resistance), but unfortunatly it doesnt work well. I would be very happy if someone can help me in this regard.
I am sharing my code to be more precise.
X0 = [0; 0; 0]; % initial conditions
r = linspace(1, 50, 100); % resistance range
freq_vector = linspace(1, 20, 20); % frequency range
options = odeset(‘RelTol’, 1e-6, ‘AbsTol’, 1e-6);
% Initialize array to store average power values over frequency range
avgpower_total = zeros(length(freq_vector), length(r));
% Loop over frequencies
for j = 1:length(freq_vector)
w = freq_vector(j); % Current frequency
avgpower = zeros(1, length(r));
tspan = linspace(0,100, 100);
F = 0.8*9.81*cos(w*tspan);
% Loop over resistance values
for i = 1:length(r)
R = r(i); % Current resistance
% Solve ODE using ode45
[t, X] = ode45(@(t, X) VoltResistFun2(t, X, R,F, tspan,w), tspan, X0, options);
% Extraction of voltage value from the solution
Voltage = X(:, 3); % Third state variable
% Calculation of power
P = Voltage.^2 / R;
% Integrated power over time
IntegratedPower = trapz(t, P);
% Calculate average power
avgpower(i) = IntegratedPower / (t(end) – t(1)); % Divide by total time
end
% Store average power for current frequency
avgpower_total(j, 🙂 = avgpower;
% keyboard
end
% Integrate the power over the frequency range and divide by the range
avg_power_over_range = trapz(freq_vector, avgpower_total, 1) / (freq_vector(end) – freq_vector(1));
function dxdt = VoltResistFun2(t, X, R,F, tspan,w)
% Parameters
g = 9.81; % m/s^2
M = 0.01; % proof mass
t_p = 0.01; % thickness
A_p = 0.0001; % area
M_p = 0.00075; % proof mass
E_33 = 1.137e-8; % permittivity
k_33 = 0.75; % electromechanical coupling
e_33 = -1; % value of e_33
m = M + (1/3) * M_p;
C_p = E_33 * A_p / t_p;
theta = -(e_33 * A_p / t_p);
alpha = 1 ./ R;
w0 = 2 * pi * 20;
% F = 0.8 * g*cos(t); % external force
z = 0.02; % damping coefficient
f = interp1(tspan, F.*cos(w*t), t);
% Define state-space matrices
A = [0 1 0; -w0^2 -2*z*w0 -theta/m; 0 theta/C_p -alpha/C_p];
B = [0; -1; 0];
% Calculate derivative of state vector
% keyboard
dxdt = A(X – B*f);
% Ensure dxdt is a column vector
% dxdt = dxdt(:);
end
Note :
eqution1 = x" +2zetaw_0x’ +w_0^2x+theta/m= -F
F is defined
1/R = alpha
Thanks in advanceI need help in solving the second order coupled differential equations as I’m unable to get the right solution after spending many days.
I want to solve it with two for loops (One for frequency and other for Resistance), but unfortunatly it doesnt work well. I would be very happy if someone can help me in this regard.
I am sharing my code to be more precise.
X0 = [0; 0; 0]; % initial conditions
r = linspace(1, 50, 100); % resistance range
freq_vector = linspace(1, 20, 20); % frequency range
options = odeset(‘RelTol’, 1e-6, ‘AbsTol’, 1e-6);
% Initialize array to store average power values over frequency range
avgpower_total = zeros(length(freq_vector), length(r));
% Loop over frequencies
for j = 1:length(freq_vector)
w = freq_vector(j); % Current frequency
avgpower = zeros(1, length(r));
tspan = linspace(0,100, 100);
F = 0.8*9.81*cos(w*tspan);
% Loop over resistance values
for i = 1:length(r)
R = r(i); % Current resistance
% Solve ODE using ode45
[t, X] = ode45(@(t, X) VoltResistFun2(t, X, R,F, tspan,w), tspan, X0, options);
% Extraction of voltage value from the solution
Voltage = X(:, 3); % Third state variable
% Calculation of power
P = Voltage.^2 / R;
% Integrated power over time
IntegratedPower = trapz(t, P);
% Calculate average power
avgpower(i) = IntegratedPower / (t(end) – t(1)); % Divide by total time
end
% Store average power for current frequency
avgpower_total(j, 🙂 = avgpower;
% keyboard
end
% Integrate the power over the frequency range and divide by the range
avg_power_over_range = trapz(freq_vector, avgpower_total, 1) / (freq_vector(end) – freq_vector(1));
function dxdt = VoltResistFun2(t, X, R,F, tspan,w)
% Parameters
g = 9.81; % m/s^2
M = 0.01; % proof mass
t_p = 0.01; % thickness
A_p = 0.0001; % area
M_p = 0.00075; % proof mass
E_33 = 1.137e-8; % permittivity
k_33 = 0.75; % electromechanical coupling
e_33 = -1; % value of e_33
m = M + (1/3) * M_p;
C_p = E_33 * A_p / t_p;
theta = -(e_33 * A_p / t_p);
alpha = 1 ./ R;
w0 = 2 * pi * 20;
% F = 0.8 * g*cos(t); % external force
z = 0.02; % damping coefficient
f = interp1(tspan, F.*cos(w*t), t);
% Define state-space matrices
A = [0 1 0; -w0^2 -2*z*w0 -theta/m; 0 theta/C_p -alpha/C_p];
B = [0; -1; 0];
% Calculate derivative of state vector
% keyboard
dxdt = A(X – B*f);
% Ensure dxdt is a column vector
% dxdt = dxdt(:);
end
Note :
eqution1 = x" +2zetaw_0x’ +w_0^2x+theta/m= -F
F is defined
1/R = alpha
Thanks in advance I need help in solving the second order coupled differential equations as I’m unable to get the right solution after spending many days.
I want to solve it with two for loops (One for frequency and other for Resistance), but unfortunatly it doesnt work well. I would be very happy if someone can help me in this regard.
I am sharing my code to be more precise.
X0 = [0; 0; 0]; % initial conditions
r = linspace(1, 50, 100); % resistance range
freq_vector = linspace(1, 20, 20); % frequency range
options = odeset(‘RelTol’, 1e-6, ‘AbsTol’, 1e-6);
% Initialize array to store average power values over frequency range
avgpower_total = zeros(length(freq_vector), length(r));
% Loop over frequencies
for j = 1:length(freq_vector)
w = freq_vector(j); % Current frequency
avgpower = zeros(1, length(r));
tspan = linspace(0,100, 100);
F = 0.8*9.81*cos(w*tspan);
% Loop over resistance values
for i = 1:length(r)
R = r(i); % Current resistance
% Solve ODE using ode45
[t, X] = ode45(@(t, X) VoltResistFun2(t, X, R,F, tspan,w), tspan, X0, options);
% Extraction of voltage value from the solution
Voltage = X(:, 3); % Third state variable
% Calculation of power
P = Voltage.^2 / R;
% Integrated power over time
IntegratedPower = trapz(t, P);
% Calculate average power
avgpower(i) = IntegratedPower / (t(end) – t(1)); % Divide by total time
end
% Store average power for current frequency
avgpower_total(j, 🙂 = avgpower;
% keyboard
end
% Integrate the power over the frequency range and divide by the range
avg_power_over_range = trapz(freq_vector, avgpower_total, 1) / (freq_vector(end) – freq_vector(1));
function dxdt = VoltResistFun2(t, X, R,F, tspan,w)
% Parameters
g = 9.81; % m/s^2
M = 0.01; % proof mass
t_p = 0.01; % thickness
A_p = 0.0001; % area
M_p = 0.00075; % proof mass
E_33 = 1.137e-8; % permittivity
k_33 = 0.75; % electromechanical coupling
e_33 = -1; % value of e_33
m = M + (1/3) * M_p;
C_p = E_33 * A_p / t_p;
theta = -(e_33 * A_p / t_p);
alpha = 1 ./ R;
w0 = 2 * pi * 20;
% F = 0.8 * g*cos(t); % external force
z = 0.02; % damping coefficient
f = interp1(tspan, F.*cos(w*t), t);
% Define state-space matrices
A = [0 1 0; -w0^2 -2*z*w0 -theta/m; 0 theta/C_p -alpha/C_p];
B = [0; -1; 0];
% Calculate derivative of state vector
% keyboard
dxdt = A(X – B*f);
% Ensure dxdt is a column vector
% dxdt = dxdt(:);
end
Note :
eqution1 = x" +2zetaw_0x’ +w_0^2x+theta/m= -F
F is defined
1/R = alpha
Thanks in advance ode45, second order coupled equations, energy harvester, avg power calculation MATLAB Answers — New Questions
Cannot load python library with R2024a: API functions not available
I updated MATLAB to R2024a a few days ago, and I’m now trying the basic example to call Python from MATLAB here: https://ch.mathworks.com/help/matlab/call-python-libraries.html
I tried the following:
>> py.list
and see this error immediately:
Python API functions are not available.
I have the PythonEnvironment setup as follows:
>> pyversion(‘/Users/…/opt/anaconda3/bin/python’)
>> pyenv
ans =
PythonEnvironment with properties:
Version: "3.9"
Executable: "/Users/…/opt/anaconda3/bin/python"
Library: "/Users/…/opt/anaconda3/lib/libpython3.9.dylib"
Home: "/Users/…/opt/anaconda3"
Status: NotLoaded
ExecutionMode: OutOfProcess
It appears I cannot get the status to turn to "Loaded", which used to happen earlier with R2023b. I am also using the apple silicon native version (maca64) of MATLAB, which could also be a factor here. I read https://ch.mathworks.com/matlabcentral/answers/1977529-how-to-use-python-from-matlab-on-mac-with-apple-silicon and have installed Amazon Corretto 11 from here: https://ch.mathworks.com/support/requirements/apple-silicon.html
Any idea how I could get the Python API to load? Many thanks in advance for your help!I updated MATLAB to R2024a a few days ago, and I’m now trying the basic example to call Python from MATLAB here: https://ch.mathworks.com/help/matlab/call-python-libraries.html
I tried the following:
>> py.list
and see this error immediately:
Python API functions are not available.
I have the PythonEnvironment setup as follows:
>> pyversion(‘/Users/…/opt/anaconda3/bin/python’)
>> pyenv
ans =
PythonEnvironment with properties:
Version: "3.9"
Executable: "/Users/…/opt/anaconda3/bin/python"
Library: "/Users/…/opt/anaconda3/lib/libpython3.9.dylib"
Home: "/Users/…/opt/anaconda3"
Status: NotLoaded
ExecutionMode: OutOfProcess
It appears I cannot get the status to turn to "Loaded", which used to happen earlier with R2023b. I am also using the apple silicon native version (maca64) of MATLAB, which could also be a factor here. I read https://ch.mathworks.com/matlabcentral/answers/1977529-how-to-use-python-from-matlab-on-mac-with-apple-silicon and have installed Amazon Corretto 11 from here: https://ch.mathworks.com/support/requirements/apple-silicon.html
Any idea how I could get the Python API to load? Many thanks in advance for your help! I updated MATLAB to R2024a a few days ago, and I’m now trying the basic example to call Python from MATLAB here: https://ch.mathworks.com/help/matlab/call-python-libraries.html
I tried the following:
>> py.list
and see this error immediately:
Python API functions are not available.
I have the PythonEnvironment setup as follows:
>> pyversion(‘/Users/…/opt/anaconda3/bin/python’)
>> pyenv
ans =
PythonEnvironment with properties:
Version: "3.9"
Executable: "/Users/…/opt/anaconda3/bin/python"
Library: "/Users/…/opt/anaconda3/lib/libpython3.9.dylib"
Home: "/Users/…/opt/anaconda3"
Status: NotLoaded
ExecutionMode: OutOfProcess
It appears I cannot get the status to turn to "Loaded", which used to happen earlier with R2023b. I am also using the apple silicon native version (maca64) of MATLAB, which could also be a factor here. I read https://ch.mathworks.com/matlabcentral/answers/1977529-how-to-use-python-from-matlab-on-mac-with-apple-silicon and have installed Amazon Corretto 11 from here: https://ch.mathworks.com/support/requirements/apple-silicon.html
Any idea how I could get the Python API to load? Many thanks in advance for your help! python MATLAB Answers — New Questions
The “parfeval” command is not working on the addPath function
I have a cell array that contains the folder paths that I want to add to the project path. Initially, I used the clasic for loop; however, as the cell array contains more than 2K folder paths, this clasic for loop takes more than 2 hours to complete the task. To overcome this issue, I have used the "parfeval" command in the for loop, but I am getting the error mentioned below.
Error Details:
Caused by:
Error using matlab.project.Project/addPath
Dot indexing is not supported for variables of this type.
How can i resolve this issue?
Below is my code:
%get current project
proj_modelBasedDesignCodeGen = currentProject;
% List of folder paths
FolderPaths = {‘componentsComponentA’, ‘componentsComponentB’, ‘componentsComponentC’, ‘componentsComponentD’};
% Preallocate cell array to store output
futures = cell(1, numel(FolderPaths));
% Loop through each folder path and submit a task to the parallel pool
for i = 1:numel(FolderPaths)
inputFolder = FolderPaths{i};
% Use parfeval to asynchronously evaluate function addPath
futures{i}= parfeval(@addPath, 1,proj_modelBasedDesignCodeGen, inputFolder); % 1 is the number of output arguments
endI have a cell array that contains the folder paths that I want to add to the project path. Initially, I used the clasic for loop; however, as the cell array contains more than 2K folder paths, this clasic for loop takes more than 2 hours to complete the task. To overcome this issue, I have used the "parfeval" command in the for loop, but I am getting the error mentioned below.
Error Details:
Caused by:
Error using matlab.project.Project/addPath
Dot indexing is not supported for variables of this type.
How can i resolve this issue?
Below is my code:
%get current project
proj_modelBasedDesignCodeGen = currentProject;
% List of folder paths
FolderPaths = {‘componentsComponentA’, ‘componentsComponentB’, ‘componentsComponentC’, ‘componentsComponentD’};
% Preallocate cell array to store output
futures = cell(1, numel(FolderPaths));
% Loop through each folder path and submit a task to the parallel pool
for i = 1:numel(FolderPaths)
inputFolder = FolderPaths{i};
% Use parfeval to asynchronously evaluate function addPath
futures{i}= parfeval(@addPath, 1,proj_modelBasedDesignCodeGen, inputFolder); % 1 is the number of output arguments
end I have a cell array that contains the folder paths that I want to add to the project path. Initially, I used the clasic for loop; however, as the cell array contains more than 2K folder paths, this clasic for loop takes more than 2 hours to complete the task. To overcome this issue, I have used the "parfeval" command in the for loop, but I am getting the error mentioned below.
Error Details:
Caused by:
Error using matlab.project.Project/addPath
Dot indexing is not supported for variables of this type.
How can i resolve this issue?
Below is my code:
%get current project
proj_modelBasedDesignCodeGen = currentProject;
% List of folder paths
FolderPaths = {‘componentsComponentA’, ‘componentsComponentB’, ‘componentsComponentC’, ‘componentsComponentD’};
% Preallocate cell array to store output
futures = cell(1, numel(FolderPaths));
% Loop through each folder path and submit a task to the parallel pool
for i = 1:numel(FolderPaths)
inputFolder = FolderPaths{i};
% Use parfeval to asynchronously evaluate function addPath
futures{i}= parfeval(@addPath, 1,proj_modelBasedDesignCodeGen, inputFolder); % 1 is the number of output arguments
end parallel computing MATLAB Answers — New Questions
How to re-generate script from already built workspace
I am having some difficulty in finding my script for several differential equations solution. Please someone help me in how to regenerate MATLAB script from workspace because i have variable workspace. when i first loaded the script it was saved in workspace format and i can not retrieve the script again because of lack of experienceI am having some difficulty in finding my script for several differential equations solution. Please someone help me in how to regenerate MATLAB script from workspace because i have variable workspace. when i first loaded the script it was saved in workspace format and i can not retrieve the script again because of lack of experience I am having some difficulty in finding my script for several differential equations solution. Please someone help me in how to regenerate MATLAB script from workspace because i have variable workspace. when i first loaded the script it was saved in workspace format and i can not retrieve the script again because of lack of experience script, workspace, homework MATLAB Answers — New Questions
Does MPC without constraints also solves QP problem? How to check QP solver within matlab?
Hi, I have one confusion about MPC controller, usually the cost function is minimized by solving QP problem at each time step.
Is MPC solves QP problem only when we have constraints or unconstraint’s MPC the optimal U is find by solving QP problem? and How to check in MATLAB that MPC is solving QP or not?Hi, I have one confusion about MPC controller, usually the cost function is minimized by solving QP problem at each time step.
Is MPC solves QP problem only when we have constraints or unconstraint’s MPC the optimal U is find by solving QP problem? and How to check in MATLAB that MPC is solving QP or not? Hi, I have one confusion about MPC controller, usually the cost function is minimized by solving QP problem at each time step.
Is MPC solves QP problem only when we have constraints or unconstraint’s MPC the optimal U is find by solving QP problem? and How to check in MATLAB that MPC is solving QP or not? model predictive control, mpc, optimization problem, qp problem MATLAB Answers — New Questions
How to get the screen size where the current uifigure is located?
I want to capture the current uifigure while running, and I have test the function with this demo code.
there is one issue should to be solved: the Y postion from the demo code is base on the left top corner, but the Y position get from uifigure.positon is based on the left bottom corner. So it need to be converted to suit for the code. So my question is how to get the screen size where the active uifigure is located in the multiple screens?
(get(0,’screensize’) can’t get the secondary screen size)
other answers which can capture the active uifigure is also welcome!I want to capture the current uifigure while running, and I have test the function with this demo code.
there is one issue should to be solved: the Y postion from the demo code is base on the left top corner, but the Y position get from uifigure.positon is based on the left bottom corner. So it need to be converted to suit for the code. So my question is how to get the screen size where the active uifigure is located in the multiple screens?
(get(0,’screensize’) can’t get the secondary screen size)
other answers which can capture the active uifigure is also welcome! I want to capture the current uifigure while running, and I have test the function with this demo code.
there is one issue should to be solved: the Y postion from the demo code is base on the left top corner, but the Y position get from uifigure.positon is based on the left bottom corner. So it need to be converted to suit for the code. So my question is how to get the screen size where the active uifigure is located in the multiple screens?
(get(0,’screensize’) can’t get the secondary screen size)
other answers which can capture the active uifigure is also welcome! capture, screen size, get screen size MATLAB Answers — New Questions
Getting screen resolution on Windows 11
Are there other functions I can use to get the monitor resolution on Windows 11? I use this code:
set(0, ‘units’, ‘pixels’); % Sets the units of your root object (screen) to pixels
Pix_SS = get(0, ‘screensize’); % Obtains the pixel information
And all I get is 1920 x 1080? This Dell U4320Q is 3840 x 2160 as shown in this Settings > Display resolution:Are there other functions I can use to get the monitor resolution on Windows 11? I use this code:
set(0, ‘units’, ‘pixels’); % Sets the units of your root object (screen) to pixels
Pix_SS = get(0, ‘screensize’); % Obtains the pixel information
And all I get is 1920 x 1080? This Dell U4320Q is 3840 x 2160 as shown in this Settings > Display resolution: Are there other functions I can use to get the monitor resolution on Windows 11? I use this code:
set(0, ‘units’, ‘pixels’); % Sets the units of your root object (screen) to pixels
Pix_SS = get(0, ‘screensize’); % Obtains the pixel information
And all I get is 1920 x 1080? This Dell U4320Q is 3840 x 2160 as shown in this Settings > Display resolution: windows 11, monitor resolution, 4k MATLAB Answers — New Questions
Problem with the fitting
% Define the data
dr_data = [2.5453, 0.042123; 5.0907, 0.075326; 7.636, 0.059506; 10.1813, 0.071553; 12.7267, 0.071365; 15.272, 0.067195; 17.8173, 0.046372; 20.3627, 0.043397; 22.908, 0.017179; 25.4533, -0.0063329; 27.9987, -0.030789; 30.544, -0.047569; 33.0893, -0.089512; 35.6347, -0.080675; 38.18, -0.089138; 40.7253, -0.1102; 43.2707, -0.12061; 45.816, -0.11857; 48.3613, -0.11955; 50.9067, -0.10803; 53.452, -0.10462; 55.9973, -0.099548; 58.5427, -0.097164; 61.088, -0.09994; 63.6333, -0.077017; 66.1787, -0.062839; 68.724, -0.048422; 71.2693, -0.03686; 73.8147, -0.01469; 76.3, 0];
dtheta_data = [2.5453, -0.099251; 5.0907, -0.16064; 7.636, -0.21858; 10.1813, -0.18965; 12.7267, -0.16996; 15.272, -0.18172; 17.8173, -0.15029; 20.3627, -0.12541; 22.908, -0.082786; 25.4533, -0.0071716; 27.9987, 0.03695; 30.544, 0.089002; 33.0893, 0.12873; 35.6347, 0.13092; 38.18, 0.13908; 40.7253, 0.17211; 43.2707, 0.16686; 45.816, 0.15826; 48.3613, 0.14872; 50.9067, 0.15295; 53.452, 0.12677; 55.9973, 0.10964; 58.5427, 0.10223; 61.088, 0.10951; 63.6333, 0.088493; 66.1787, 0.068903; 68.724, 0.054396; 71.2693, 0.035731; 73.8147, 0.030172; 76.3, 0];
% Define initial guess values for parameters
lambda_init = 2.9;
kappa_init = 0.066;
theta_k_init = 0.45;
R_init = 2400;
rout = 76.3;
% Define c
c = sqrt(4 – (rout/R_init)^2);
% Define initial conditions as parameters
a0_bound = 0.0;
b0_bound = 0.0;
c0_bound = 0.0;
d0_bound = 0.0;
% Objective function
objective_function = @(params) compute_error(params, dr_data, dtheta_data, R_init, rout, c);
% Optimization options
options = optimoptions(@fminunc,’Display’,’iter’,’Algorithm’,’quasi-newton’,’MaxFunctionEvaluations’,10000,’MaxIterations’,10000);
% Perform optimization
best_params = fminunc(objective_function, [lambda_init, kappa_init, theta_k_init, a0_bound, b0_bound, c0_bound, d0_bound], options);
% Display best parameters
disp("Optimal Parameters:");
disp("Lambda: " + best_params(1));
disp("Kappa: " + best_params(2));
disp("Theta_k: " + best_params(3));
% Plot the results with the best parameters
[~, dr_fit, dtheta_fit] = compute_error(best_params, dr_data, dtheta_data, R_init, rout, c);
figure;
subplot(2,1,1);
plot(dr_data(:,1), dr_data(:,2), ‘ro’, dr_data(:,1), dr_fit, ‘b-‘);
xlabel(‘r’);
ylabel(‘dr(r)’);
title(‘Fit of dr(r) vs r’);
legend(‘Data’, ‘Fit’);
subplot(2,1,2);
plot(dtheta_data(:,1), dtheta_data(:,2), ‘ro’, dtheta_data(:,1), dtheta_fit, ‘b-‘);
xlabel(‘r’);
ylabel(‘dtheta(r)’);
title(‘Fit of dtheta(r) vs r’);
legend(‘Data’, ‘Fit’);
% Error computation function
function [error, dr_fit, dtheta_fit] = compute_error(params, dr_data, dtheta_data, R, rout, c)
lambda = params(1);
kappa = params(2);
theta_k = params(3);
a0 = params(4);
b0 = params(5);
c0 = params(6);
d0 = params(7);
% Define the function for the differential equations
f = @(r, y) [y(2); ((-1*(lambda+1)*(r*y(2)+y(1)))+(1/r)*((kappa*r^2*cos(theta_k)-(lambda+1))*y(1))-(kappa*r*y(3)*sin(theta_k))+(-16*lambda*r^2)/(c^4*R^2)); y(4); ((-1*(r*y(4)+y(3)))+(1/r)*((kappa*r^2*cos(theta_k)-1)*y(3))+(kappa*r*y(1)*sin(theta_k)))];
% Solve the differential equations using ode45
%r_span = [min(dr_data(:,1)), max(dr_data(:,1))]; % Span of r values from the data
r_span = dr_data(:,1);
[~, sol] = ode45(f, r_span, [a0, b0, c0, d0]);
% Extract solutions
dr_fit = sol(:,1);
dtheta_fit = sol(:,3);
% Compute error (sum of squared differences)
error_dr = sum((dr_data(:,2) – dr_fit).^2);
error_dtheta = sum((dtheta_data(:,2) – dtheta_fit).^2);
% Total error
error = error_dr + error_dtheta;
end
I want to fit the simulated d_r(r) vs r and d_theta(r) vs r values with the above mentioned coupled differential eqns by usuing the fitting parameters lamda, kappa, theta_k. For the sake of good fitting one can use boundary conditions as parameter also. However getting errors. Please help me to solve those errors and fitting those data.
The mathematical eqns and details are below:
$kappa>0$ and $0 leq theta_{k}<frac{pi}{2}$ describe the screening. The boundary conditions are $vec{d}(0)=vec{d}left(r_{text {out }}right)=(0,0)$.
For
$$
c=sqrt{4-left(frac{r_{text {out }}}{R}right)^{2}}
$$
Furthermore,
begin{equation}
(lambda+1)left(r d_{r}^{prime prime}(r)+d_{r}^{prime}(r)right) &+ frac{1}{r}left(kappa r^{2} cos theta_{k}-(lambda+1)right) d_{r}(r)\
& – kappa r d_{theta}(r) sin theta_{k} = -frac{16 lambda r^{2}}{c^{4} R^{2}}
end{equation}
begin{equation}
r d_{theta}^{prime prime}(r) + d_{theta}^{prime}(r) &+ frac{1}{r}left(kappa r^{2} cos theta_{k}-1right) d_{theta}(r)\
& + kappa r d_{r}(r) sin theta_{k} = 0
end{equation}% Define the data
dr_data = [2.5453, 0.042123; 5.0907, 0.075326; 7.636, 0.059506; 10.1813, 0.071553; 12.7267, 0.071365; 15.272, 0.067195; 17.8173, 0.046372; 20.3627, 0.043397; 22.908, 0.017179; 25.4533, -0.0063329; 27.9987, -0.030789; 30.544, -0.047569; 33.0893, -0.089512; 35.6347, -0.080675; 38.18, -0.089138; 40.7253, -0.1102; 43.2707, -0.12061; 45.816, -0.11857; 48.3613, -0.11955; 50.9067, -0.10803; 53.452, -0.10462; 55.9973, -0.099548; 58.5427, -0.097164; 61.088, -0.09994; 63.6333, -0.077017; 66.1787, -0.062839; 68.724, -0.048422; 71.2693, -0.03686; 73.8147, -0.01469; 76.3, 0];
dtheta_data = [2.5453, -0.099251; 5.0907, -0.16064; 7.636, -0.21858; 10.1813, -0.18965; 12.7267, -0.16996; 15.272, -0.18172; 17.8173, -0.15029; 20.3627, -0.12541; 22.908, -0.082786; 25.4533, -0.0071716; 27.9987, 0.03695; 30.544, 0.089002; 33.0893, 0.12873; 35.6347, 0.13092; 38.18, 0.13908; 40.7253, 0.17211; 43.2707, 0.16686; 45.816, 0.15826; 48.3613, 0.14872; 50.9067, 0.15295; 53.452, 0.12677; 55.9973, 0.10964; 58.5427, 0.10223; 61.088, 0.10951; 63.6333, 0.088493; 66.1787, 0.068903; 68.724, 0.054396; 71.2693, 0.035731; 73.8147, 0.030172; 76.3, 0];
% Define initial guess values for parameters
lambda_init = 2.9;
kappa_init = 0.066;
theta_k_init = 0.45;
R_init = 2400;
rout = 76.3;
% Define c
c = sqrt(4 – (rout/R_init)^2);
% Define initial conditions as parameters
a0_bound = 0.0;
b0_bound = 0.0;
c0_bound = 0.0;
d0_bound = 0.0;
% Objective function
objective_function = @(params) compute_error(params, dr_data, dtheta_data, R_init, rout, c);
% Optimization options
options = optimoptions(@fminunc,’Display’,’iter’,’Algorithm’,’quasi-newton’,’MaxFunctionEvaluations’,10000,’MaxIterations’,10000);
% Perform optimization
best_params = fminunc(objective_function, [lambda_init, kappa_init, theta_k_init, a0_bound, b0_bound, c0_bound, d0_bound], options);
% Display best parameters
disp("Optimal Parameters:");
disp("Lambda: " + best_params(1));
disp("Kappa: " + best_params(2));
disp("Theta_k: " + best_params(3));
% Plot the results with the best parameters
[~, dr_fit, dtheta_fit] = compute_error(best_params, dr_data, dtheta_data, R_init, rout, c);
figure;
subplot(2,1,1);
plot(dr_data(:,1), dr_data(:,2), ‘ro’, dr_data(:,1), dr_fit, ‘b-‘);
xlabel(‘r’);
ylabel(‘dr(r)’);
title(‘Fit of dr(r) vs r’);
legend(‘Data’, ‘Fit’);
subplot(2,1,2);
plot(dtheta_data(:,1), dtheta_data(:,2), ‘ro’, dtheta_data(:,1), dtheta_fit, ‘b-‘);
xlabel(‘r’);
ylabel(‘dtheta(r)’);
title(‘Fit of dtheta(r) vs r’);
legend(‘Data’, ‘Fit’);
% Error computation function
function [error, dr_fit, dtheta_fit] = compute_error(params, dr_data, dtheta_data, R, rout, c)
lambda = params(1);
kappa = params(2);
theta_k = params(3);
a0 = params(4);
b0 = params(5);
c0 = params(6);
d0 = params(7);
% Define the function for the differential equations
f = @(r, y) [y(2); ((-1*(lambda+1)*(r*y(2)+y(1)))+(1/r)*((kappa*r^2*cos(theta_k)-(lambda+1))*y(1))-(kappa*r*y(3)*sin(theta_k))+(-16*lambda*r^2)/(c^4*R^2)); y(4); ((-1*(r*y(4)+y(3)))+(1/r)*((kappa*r^2*cos(theta_k)-1)*y(3))+(kappa*r*y(1)*sin(theta_k)))];
% Solve the differential equations using ode45
%r_span = [min(dr_data(:,1)), max(dr_data(:,1))]; % Span of r values from the data
r_span = dr_data(:,1);
[~, sol] = ode45(f, r_span, [a0, b0, c0, d0]);
% Extract solutions
dr_fit = sol(:,1);
dtheta_fit = sol(:,3);
% Compute error (sum of squared differences)
error_dr = sum((dr_data(:,2) – dr_fit).^2);
error_dtheta = sum((dtheta_data(:,2) – dtheta_fit).^2);
% Total error
error = error_dr + error_dtheta;
end
I want to fit the simulated d_r(r) vs r and d_theta(r) vs r values with the above mentioned coupled differential eqns by usuing the fitting parameters lamda, kappa, theta_k. For the sake of good fitting one can use boundary conditions as parameter also. However getting errors. Please help me to solve those errors and fitting those data.
The mathematical eqns and details are below:
$kappa>0$ and $0 leq theta_{k}<frac{pi}{2}$ describe the screening. The boundary conditions are $vec{d}(0)=vec{d}left(r_{text {out }}right)=(0,0)$.
For
$$
c=sqrt{4-left(frac{r_{text {out }}}{R}right)^{2}}
$$
Furthermore,
begin{equation}
(lambda+1)left(r d_{r}^{prime prime}(r)+d_{r}^{prime}(r)right) &+ frac{1}{r}left(kappa r^{2} cos theta_{k}-(lambda+1)right) d_{r}(r)\
& – kappa r d_{theta}(r) sin theta_{k} = -frac{16 lambda r^{2}}{c^{4} R^{2}}
end{equation}
begin{equation}
r d_{theta}^{prime prime}(r) + d_{theta}^{prime}(r) &+ frac{1}{r}left(kappa r^{2} cos theta_{k}-1right) d_{theta}(r)\
& + kappa r d_{r}(r) sin theta_{k} = 0
end{equation} % Define the data
dr_data = [2.5453, 0.042123; 5.0907, 0.075326; 7.636, 0.059506; 10.1813, 0.071553; 12.7267, 0.071365; 15.272, 0.067195; 17.8173, 0.046372; 20.3627, 0.043397; 22.908, 0.017179; 25.4533, -0.0063329; 27.9987, -0.030789; 30.544, -0.047569; 33.0893, -0.089512; 35.6347, -0.080675; 38.18, -0.089138; 40.7253, -0.1102; 43.2707, -0.12061; 45.816, -0.11857; 48.3613, -0.11955; 50.9067, -0.10803; 53.452, -0.10462; 55.9973, -0.099548; 58.5427, -0.097164; 61.088, -0.09994; 63.6333, -0.077017; 66.1787, -0.062839; 68.724, -0.048422; 71.2693, -0.03686; 73.8147, -0.01469; 76.3, 0];
dtheta_data = [2.5453, -0.099251; 5.0907, -0.16064; 7.636, -0.21858; 10.1813, -0.18965; 12.7267, -0.16996; 15.272, -0.18172; 17.8173, -0.15029; 20.3627, -0.12541; 22.908, -0.082786; 25.4533, -0.0071716; 27.9987, 0.03695; 30.544, 0.089002; 33.0893, 0.12873; 35.6347, 0.13092; 38.18, 0.13908; 40.7253, 0.17211; 43.2707, 0.16686; 45.816, 0.15826; 48.3613, 0.14872; 50.9067, 0.15295; 53.452, 0.12677; 55.9973, 0.10964; 58.5427, 0.10223; 61.088, 0.10951; 63.6333, 0.088493; 66.1787, 0.068903; 68.724, 0.054396; 71.2693, 0.035731; 73.8147, 0.030172; 76.3, 0];
% Define initial guess values for parameters
lambda_init = 2.9;
kappa_init = 0.066;
theta_k_init = 0.45;
R_init = 2400;
rout = 76.3;
% Define c
c = sqrt(4 – (rout/R_init)^2);
% Define initial conditions as parameters
a0_bound = 0.0;
b0_bound = 0.0;
c0_bound = 0.0;
d0_bound = 0.0;
% Objective function
objective_function = @(params) compute_error(params, dr_data, dtheta_data, R_init, rout, c);
% Optimization options
options = optimoptions(@fminunc,’Display’,’iter’,’Algorithm’,’quasi-newton’,’MaxFunctionEvaluations’,10000,’MaxIterations’,10000);
% Perform optimization
best_params = fminunc(objective_function, [lambda_init, kappa_init, theta_k_init, a0_bound, b0_bound, c0_bound, d0_bound], options);
% Display best parameters
disp("Optimal Parameters:");
disp("Lambda: " + best_params(1));
disp("Kappa: " + best_params(2));
disp("Theta_k: " + best_params(3));
% Plot the results with the best parameters
[~, dr_fit, dtheta_fit] = compute_error(best_params, dr_data, dtheta_data, R_init, rout, c);
figure;
subplot(2,1,1);
plot(dr_data(:,1), dr_data(:,2), ‘ro’, dr_data(:,1), dr_fit, ‘b-‘);
xlabel(‘r’);
ylabel(‘dr(r)’);
title(‘Fit of dr(r) vs r’);
legend(‘Data’, ‘Fit’);
subplot(2,1,2);
plot(dtheta_data(:,1), dtheta_data(:,2), ‘ro’, dtheta_data(:,1), dtheta_fit, ‘b-‘);
xlabel(‘r’);
ylabel(‘dtheta(r)’);
title(‘Fit of dtheta(r) vs r’);
legend(‘Data’, ‘Fit’);
% Error computation function
function [error, dr_fit, dtheta_fit] = compute_error(params, dr_data, dtheta_data, R, rout, c)
lambda = params(1);
kappa = params(2);
theta_k = params(3);
a0 = params(4);
b0 = params(5);
c0 = params(6);
d0 = params(7);
% Define the function for the differential equations
f = @(r, y) [y(2); ((-1*(lambda+1)*(r*y(2)+y(1)))+(1/r)*((kappa*r^2*cos(theta_k)-(lambda+1))*y(1))-(kappa*r*y(3)*sin(theta_k))+(-16*lambda*r^2)/(c^4*R^2)); y(4); ((-1*(r*y(4)+y(3)))+(1/r)*((kappa*r^2*cos(theta_k)-1)*y(3))+(kappa*r*y(1)*sin(theta_k)))];
% Solve the differential equations using ode45
%r_span = [min(dr_data(:,1)), max(dr_data(:,1))]; % Span of r values from the data
r_span = dr_data(:,1);
[~, sol] = ode45(f, r_span, [a0, b0, c0, d0]);
% Extract solutions
dr_fit = sol(:,1);
dtheta_fit = sol(:,3);
% Compute error (sum of squared differences)
error_dr = sum((dr_data(:,2) – dr_fit).^2);
error_dtheta = sum((dtheta_data(:,2) – dtheta_fit).^2);
% Total error
error = error_dr + error_dtheta;
end
I want to fit the simulated d_r(r) vs r and d_theta(r) vs r values with the above mentioned coupled differential eqns by usuing the fitting parameters lamda, kappa, theta_k. For the sake of good fitting one can use boundary conditions as parameter also. However getting errors. Please help me to solve those errors and fitting those data.
The mathematical eqns and details are below:
$kappa>0$ and $0 leq theta_{k}<frac{pi}{2}$ describe the screening. The boundary conditions are $vec{d}(0)=vec{d}left(r_{text {out }}right)=(0,0)$.
For
$$
c=sqrt{4-left(frac{r_{text {out }}}{R}right)^{2}}
$$
Furthermore,
begin{equation}
(lambda+1)left(r d_{r}^{prime prime}(r)+d_{r}^{prime}(r)right) &+ frac{1}{r}left(kappa r^{2} cos theta_{k}-(lambda+1)right) d_{r}(r)\
& – kappa r d_{theta}(r) sin theta_{k} = -frac{16 lambda r^{2}}{c^{4} R^{2}}
end{equation}
begin{equation}
r d_{theta}^{prime prime}(r) + d_{theta}^{prime}(r) &+ frac{1}{r}left(kappa r^{2} cos theta_{k}-1right) d_{theta}(r)\
& + kappa r d_{r}(r) sin theta_{k} = 0
end{equation} curve fitting, parametes, solve eqns MATLAB Answers — New Questions