script no longer works after update to 2023b (related to ode45)
17 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Patrick D Sinko
am 25 Okt. 2023
Kommentiert: Walter Roberson
am 26 Okt. 2023
I am solving the blasius boundary layer equations for concentration and momentum boundary layer with this script. it seems like something changed in the 'ode45' function from the update to 2023a to 2023b, because MATLAB seems to want want to use 'odeset' instead of 'odeget' but the error is occuring in the 'ode45' function. Are there any suggestions on how to fix this? Contents of this help request 1) function file for the ODE to be solved 2) loop to solve multiple initial conditions 3) Error message from MATLAB.
1) function file:
%% Defining the function for the blasius solution
% define initial conditions
function df = Blas_Sc(eta,f,Sc)
global Sc;
df = [f(2); f(3); -0.5*f(1)*f(3);f(5);-0.5*Sc*f(1)*f(5)];
end
2) loop to calculate compute boundary layers:
%% Blasius Boundary Layer Solution and CBL Solution for Flat Plate
%% Patrick D. Sinko Uppsala University, February 2, 2023
%purpose To solve the system of ODE for the hydrodynamic and concentration
%equations for flow over a flat plate
%Bisection Method Solver code modified from: Mohamad Aslani (2023).
% Compressible Similarity Solution for Flat Plate BL
% (https://www.mathworks.com/matlabcentral/fileexchange/73587-compressible
% -similarity-solution-for-flat-plate-bl), MATLAB Central File Exchange.
% Retrieved February 5, 2023.
%%
clear all
clc
close all
%% Global initialization
nu = 0.000000696;% [m^2/s] Water kinematic Viscosity at 37C
U_inf = 0.5; %m/s
global Sc;
Sc = 1; % starting value
% Initializing loop variables
step=1;
end_step=2;
% vector defining the span and interval of Schimdt number to be investigated.
SC = (Sc:step:end_step);
ff = zeros(length(SC),1);
gg = zeros(length(SC),1);
eta_f = zeros(length(SC),1);
eta_g = zeros(length(SC),1);
count_j =0; %storage variable counter
while(Sc <= end_step)
%% Initial guesses for f''(0): need not to be changed for most cases...
format long;
af = -0.0;
bf = 2.0;
ag = -20.0;
bg = 20.0;
BC_f = 1; % goal for the f' boundary condition i.e. (u/Uinf) = 1
BC_g = 1; % goal for the g boundary condition i.e. (c/cinf) = 1
%% Solver specific parameters
eta_max_ode = 10; % length of eta distance
% Maximum iterations
max_iter = 1e2;
% Iteration counter
count = 0; %convergence iteration while counter
% Tolerance for boundary condition
epsilon = 1e-6;
% Value greater than epsilon to start it off
error = 1e15;
ERR = zeros(max_iter,1); %initializing here clears the previous loop's data
ERR_f = zeros(max_iter,1);%initializing here clears the previous loop's data
ERR_g = zeros(max_iter,1);%initializing here clears the previous loop's data
%% Solve the ODE using ode45
tic
while(abs(error) > epsilon && count < max_iter)
sigmaf = (af+bf)/2; % initial guess
sigmag = (ag+bg)/2;
% sigmaf = 0.332;
f0 = [0, 0, sigmaf,0,sigmag];
options = odeset('RelTol',1e-9,'AbsTol',1e-9);
% the Hydrodynamic and concentration set of 1st order ODE
[eta,f] = ode45(@Blas_Sc,[0, eta_max_ode],f0,Sc,options);
% Perform bisection method
errorf = f(end,2) - BC_f; % error in the velocity profile
errorg = f(end,4) - BC_g; % error in the concentration profile
% error = abs(errorf) + abs(errorg) + abs(errorh)
error = abs(errorf) + abs(errorg);
% error = abs(errorg);
if(errorf > 0)
bf = sigmaf;
else
af = sigmaf;
end
if(errorg > 0)
bg = sigmag;
else
ag = sigmag;
end
ERR(count+1) = error; %tracking the convergence of the total error
ERR_f(count+1) = errorf; %tracking the convergence of the hydrodynamic error
ERR_g(count+1) = errorg; %tracking the convergence of the concentration error
% Increment the counter
count = count + 1;
end
toc
count=0; %reset the counter for next loop (next input of Sc)
count_j = count_j + 1; %pass the first calculation into the first row of the storage variables
3) The error message:
Warning: The value of local variables may have been changed to match the globals. Future versions of MATLAB will require
that you declare a variable to be global before you use that variable.
> In Blas_Sc (line 4)
In odearguments (line 92)
In ode45 (line 104)
In ABLCBL_3 (line 75)
Error using odeget
First argument must be an options structure created with ODESET.
Error in odearguments (line 124)
rtol = odeget(options,'RelTol',1e-3);
Error in ode45 (line 104)
odearguments(odeIsFuncHandle,odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);
Error in ABLCBL_3 (line 75)
[eta,f] = ode45(@Blas_Sc,[0, eta_max_ode],f0,Sc,options);
0 Kommentare
Akzeptierte Antwort
Torsten
am 25 Okt. 2023
Bearbeitet: Torsten
am 25 Okt. 2023
Why do you use ode45 for a boundary value problem ? Use bvp4c instead.
The below code doesn't give errors or warnings - at least as long as it ran under MATLAB online here.
%2) loop to calculate compute boundary layers:
%% Blasius Boundary Layer Solution and CBL Solution for Flat Plate
%% Patrick D. Sinko Uppsala University, February 2, 2023
%purpose To solve the system of ODE for the hydrodynamic and concentration
%equations for flow over a flat plate
%Bisection Method Solver code modified from: Mohamad Aslani (2023).
% Compressible Similarity Solution for Flat Plate BL
% (https://www.mathworks.com/matlabcentral/fileexchange/73587-compressible
% -similarity-solution-for-flat-plate-bl), MATLAB Central File Exchange.
% Retrieved February 5, 2023.
%%
clear all
clc
close all
%% Global initialization
nu = 0.000000696;% [m^2/s] Water kinematic Viscosity at 37C
U_inf = 0.5; %m/s
Sc = 1; % starting value
% Initializing loop variables
step=1;
end_step=2;
% vector defining the span and interval of Schimdt number to be investigated.
SC = (Sc:step:end_step);
ff = zeros(length(SC),1);
gg = zeros(length(SC),1);
eta_f = zeros(length(SC),1);
eta_g = zeros(length(SC),1);
count_j =0; %storage variable counter
while(Sc <= end_step)
%% Initial guesses for f''(0): need not to be changed for most cases...
format long;
af = -0.0;
bf = 2.0;
ag = -20.0;
bg = 20.0;
BC_f = 1; % goal for the f' boundary condition i.e. (u/Uinf) = 1
BC_g = 1; % goal for the g boundary condition i.e. (c/cinf) = 1
%% Solver specific parameters
eta_max_ode = 10; % length of eta distance
% Maximum iterations
max_iter = 1e2;
% Iteration counter
count = 0; %convergence iteration while counter
% Tolerance for boundary condition
epsilon = 1e-6;
% Value greater than epsilon to start it off
error = 1e15;
ERR = zeros(max_iter,1); %initializing here clears the previous loop's data
ERR_f = zeros(max_iter,1);%initializing here clears the previous loop's data
ERR_g = zeros(max_iter,1);%initializing here clears the previous loop's data
%% Solve the ODE using ode45
tic
while(abs(error) > epsilon && count < max_iter)
sigmaf = (af+bf)/2; % initial guess
sigmag = (ag+bg)/2;
% sigmaf = 0.332;
f0 = [0, 0, sigmaf,0,sigmag];
options = odeset('RelTol',1e-9,'AbsTol',1e-9);
% the Hydrodynamic and concentration set of 1st order ODE
[eta,f] = ode45(@(t,y)Blas_Sc(t,y,Sc),[0, eta_max_ode],f0,options);
% Perform bisection method
errorf = f(end,2) - BC_f; % error in the velocity profile
errorg = f(end,4) - BC_g; % error in the concentration profile
% error = abs(errorf) + abs(errorg) + abs(errorh)
error = abs(errorf) + abs(errorg);
% error = abs(errorg);
if(errorf > 0)
bf = sigmaf;
else
af = sigmaf;
end
if(errorg > 0)
bg = sigmag;
else
ag = sigmag;
end
ERR(count+1) = error; %tracking the convergence of the total error
ERR_f(count+1) = errorf; %tracking the convergence of the hydrodynamic error
ERR_g(count+1) = errorg; %tracking the convergence of the concentration error
% Increment the counter
count = count + 1;
end
toc
count=0; %reset the counter for next loop (next input of Sc)
count_j = count_j + 1; %pass the first calculation into the first row of the storage variables
end
%% Defining the function for the blasius solution
% define initial conditions
function df = Blas_Sc(eta,f,Sc)
df = [f(2); f(3); -0.5*f(1)*f(3);f(5);-0.5*Sc*f(1)*f(5)];
end
4 Kommentare
Weitere Antworten (1)
Walter Roberson
am 25 Okt. 2023
function df = Blas_Sc(eta,f,Sc)
global Sc;
It is an error to have a parameter with the same name as a global variable
1 Kommentar
Walter Roberson
am 26 Okt. 2023
global Sc;
Sc = 1; % starting value
That is fine. Sc is not given a value until after the global is executed.
[eta,f] = ode45(@Blas_Sc,[0, eta_max_ode],f0,Sc,options);
That undocumented syntax asks ode45 to treat Sc as the ode options, and to pass the content of options as an extra parameter to some of the calls. You should not use undocumented syntax. You should use http://www.mathworks.com/help/matlab/math/parameterizing-functions.html
function df = Blas_Sc(eta,f,Sc)
That line tells Blas_Sc to expect an optional third parameter, and to associate that third parameter with a variable named Sc
global Sc;
In current versions, that tells MATLAB to replace Sc (the parameter) with the content of the global variable Sc, unless no global variable Sc had been initialized, in which case a global variable named Sc would be created and would be given the value that had been passed into Sc, unless nothing was pased in the slot for Sc in which case the global variable Sc would be initialized as empty.
That is crazy enough to warrant a warning.
Just don't do it. If you are going to have a global variable, do not pass the same variable as a parameter. If the function does not need to write to the global, then do not use a global variable, just parameterize the function call if needed.
Error using odeget
First argument must be an options structure created with ODESET.
Yup. Remember what I said above that your ode45 call is passing Sc in the slot that ode45 expects options. The internal processing of the options parameter is going to do odeget(OPTIONS_PARAMETER) and is going to find that the value 1 has been passed as the options parameter.
If you are going to use undocumented syntax for functions, then get it right.
Siehe auch
Kategorien
Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!