Can I use a struct in an anonymous function?
I am solving a PDE via finite volume methods. As a result, I am left with a set of ODEs which I solve with ode15s. I have a bunch of constants which I use in the function defining the derivative, and I insert these by the inclusion of a struct to the function. When I run the ode15s I obtain an error which pertains to the inclusion of the struct as an input to the function. I don’t really want to include the constants within the function as that seems "messy", so I would prefer to do it via a struct.
%This is a code to solve the isothermal sintering equations.
%%—Physical parameters—
S=struct;
S.g = 9.81; %Acceleration due to gravity.
S.K = 1; %This is the Laplace constant from the sintering potential.
S.eta_0 = 1; %The shear stress of the skeleton
S.T = 3600; %This is the time of sintering.
S.M = 0.1; %This is the total mass of the rod being sintered .
S.rho_m = 1; %This is the density of the individual metallic particles.
S.N = 251; %This is the number of cells-1
L = 1; %This is the initial length of the rod
%%—Initial conditions
%Length
X=linspace(0,L,S.N); %This is the partition of the initial length
%Density
rho_0=0.5+0.1*sin(6*pi*X);
%Amount of mass
S.h=cumtrapz(X,rho_0);
%Initial velocity
u0=zeros(S.N-1,1);
% Define the system of equations as a function
y0 = [1./rho_0′; u0];
% Finite Volume Method solution loop
tspan = [0 3]; % You can adjust this based on the time range you’re interested in
[t, y] = ode15s(@ode_system, tspan, y0);
function dydt = ode_system(t, y, S)
S.M
U=S.M/(S.rho_m*S.T);
a_1=S.T/(U*S.M);
a_2=S.T*S.rho_m/S.M^2;
a_3=S.g*S.T/U;
dh=diff(S.h);
y=NaN(1,2*S.N-2);
nu = y(1:S.N-1); % u(t)
u = y(S.N:2*S.N-2); % v(t)
%Compute the left and right specific volumes:
nu_left=15*nu(1)/8-5*nu(2)/4+3*nu(3)/8;
nu_right=15*nu(S.N-3)/8-5*nu(S.N-2)/5+3*nu(S.N-1)/8;
%Compute the fluxes on the interfaces
nu_half=NaN(1:S.N);
nu_half(2:N-1)=0.5*(nu(2:S.N-1)+nu(1:S.N-2));
nu_half(1)=nu_left;
nu_half(N)=nu_right;
%Fluxes for nu equation
nu_flux=zeros(S.N,1);
nu_flux(2:S.N-1)=0.5*(u(1:S.N-2)+u(2:S.N-1));
du_dh_left=-P_L(S,nu_left)*nu_left/zeta(nu_left);
du_dh_right=-P_L(S,nu_right)*nu_right/zeta(nu_right);
nu_flux(1)=-3*du_dh_left*dh(1)/8+9*u(1)/8-u(2)/8;
nu_flux(N)=3*du_dh_right*dh(N-1)/8-9*u(N-1)/8+u(N-2)/8;
%Fluxes for u equation
u_grad=2*(u(2:S.N-2)-u(1:S.N-2))./(dh(2:S.N-1)+dh(1:S.N-2));
u_flux=NaN(S.N,1);
u_flux(2:S.N-1)=a_1*P_L(S,nu_half(2:S.N-1))+(a_2*zeta(nu_half(2:S.N-1))./nu_half(2:S.N-1)).*u_grad;
% The system of equationsY
dnu_dt = nu_flux(2:S.N)-nu_flux(1:S.N-1);
du_dt =u_flux(2:S.N)-u_flux(1:S.N-1);
% Return the derivatives
dydt = [dnu_dt; du_dt];
end
function y=P_L(S,nu)
y=S.K./nu;
end
function y=zeta(nu)
y=(2/3)*(nu).^(-2)./(nu-1);
endI am solving a PDE via finite volume methods. As a result, I am left with a set of ODEs which I solve with ode15s. I have a bunch of constants which I use in the function defining the derivative, and I insert these by the inclusion of a struct to the function. When I run the ode15s I obtain an error which pertains to the inclusion of the struct as an input to the function. I don’t really want to include the constants within the function as that seems "messy", so I would prefer to do it via a struct.
%This is a code to solve the isothermal sintering equations.
%%—Physical parameters—
S=struct;
S.g = 9.81; %Acceleration due to gravity.
S.K = 1; %This is the Laplace constant from the sintering potential.
S.eta_0 = 1; %The shear stress of the skeleton
S.T = 3600; %This is the time of sintering.
S.M = 0.1; %This is the total mass of the rod being sintered .
S.rho_m = 1; %This is the density of the individual metallic particles.
S.N = 251; %This is the number of cells-1
L = 1; %This is the initial length of the rod
%%—Initial conditions
%Length
X=linspace(0,L,S.N); %This is the partition of the initial length
%Density
rho_0=0.5+0.1*sin(6*pi*X);
%Amount of mass
S.h=cumtrapz(X,rho_0);
%Initial velocity
u0=zeros(S.N-1,1);
% Define the system of equations as a function
y0 = [1./rho_0′; u0];
% Finite Volume Method solution loop
tspan = [0 3]; % You can adjust this based on the time range you’re interested in
[t, y] = ode15s(@ode_system, tspan, y0);
function dydt = ode_system(t, y, S)
S.M
U=S.M/(S.rho_m*S.T);
a_1=S.T/(U*S.M);
a_2=S.T*S.rho_m/S.M^2;
a_3=S.g*S.T/U;
dh=diff(S.h);
y=NaN(1,2*S.N-2);
nu = y(1:S.N-1); % u(t)
u = y(S.N:2*S.N-2); % v(t)
%Compute the left and right specific volumes:
nu_left=15*nu(1)/8-5*nu(2)/4+3*nu(3)/8;
nu_right=15*nu(S.N-3)/8-5*nu(S.N-2)/5+3*nu(S.N-1)/8;
%Compute the fluxes on the interfaces
nu_half=NaN(1:S.N);
nu_half(2:N-1)=0.5*(nu(2:S.N-1)+nu(1:S.N-2));
nu_half(1)=nu_left;
nu_half(N)=nu_right;
%Fluxes for nu equation
nu_flux=zeros(S.N,1);
nu_flux(2:S.N-1)=0.5*(u(1:S.N-2)+u(2:S.N-1));
du_dh_left=-P_L(S,nu_left)*nu_left/zeta(nu_left);
du_dh_right=-P_L(S,nu_right)*nu_right/zeta(nu_right);
nu_flux(1)=-3*du_dh_left*dh(1)/8+9*u(1)/8-u(2)/8;
nu_flux(N)=3*du_dh_right*dh(N-1)/8-9*u(N-1)/8+u(N-2)/8;
%Fluxes for u equation
u_grad=2*(u(2:S.N-2)-u(1:S.N-2))./(dh(2:S.N-1)+dh(1:S.N-2));
u_flux=NaN(S.N,1);
u_flux(2:S.N-1)=a_1*P_L(S,nu_half(2:S.N-1))+(a_2*zeta(nu_half(2:S.N-1))./nu_half(2:S.N-1)).*u_grad;
% The system of equationsY
dnu_dt = nu_flux(2:S.N)-nu_flux(1:S.N-1);
du_dt =u_flux(2:S.N)-u_flux(1:S.N-1);
% Return the derivatives
dydt = [dnu_dt; du_dt];
end
function y=P_L(S,nu)
y=S.K./nu;
end
function y=zeta(nu)
y=(2/3)*(nu).^(-2)./(nu-1);
end I am solving a PDE via finite volume methods. As a result, I am left with a set of ODEs which I solve with ode15s. I have a bunch of constants which I use in the function defining the derivative, and I insert these by the inclusion of a struct to the function. When I run the ode15s I obtain an error which pertains to the inclusion of the struct as an input to the function. I don’t really want to include the constants within the function as that seems "messy", so I would prefer to do it via a struct.
%This is a code to solve the isothermal sintering equations.
%%—Physical parameters—
S=struct;
S.g = 9.81; %Acceleration due to gravity.
S.K = 1; %This is the Laplace constant from the sintering potential.
S.eta_0 = 1; %The shear stress of the skeleton
S.T = 3600; %This is the time of sintering.
S.M = 0.1; %This is the total mass of the rod being sintered .
S.rho_m = 1; %This is the density of the individual metallic particles.
S.N = 251; %This is the number of cells-1
L = 1; %This is the initial length of the rod
%%—Initial conditions
%Length
X=linspace(0,L,S.N); %This is the partition of the initial length
%Density
rho_0=0.5+0.1*sin(6*pi*X);
%Amount of mass
S.h=cumtrapz(X,rho_0);
%Initial velocity
u0=zeros(S.N-1,1);
% Define the system of equations as a function
y0 = [1./rho_0′; u0];
% Finite Volume Method solution loop
tspan = [0 3]; % You can adjust this based on the time range you’re interested in
[t, y] = ode15s(@ode_system, tspan, y0);
function dydt = ode_system(t, y, S)
S.M
U=S.M/(S.rho_m*S.T);
a_1=S.T/(U*S.M);
a_2=S.T*S.rho_m/S.M^2;
a_3=S.g*S.T/U;
dh=diff(S.h);
y=NaN(1,2*S.N-2);
nu = y(1:S.N-1); % u(t)
u = y(S.N:2*S.N-2); % v(t)
%Compute the left and right specific volumes:
nu_left=15*nu(1)/8-5*nu(2)/4+3*nu(3)/8;
nu_right=15*nu(S.N-3)/8-5*nu(S.N-2)/5+3*nu(S.N-1)/8;
%Compute the fluxes on the interfaces
nu_half=NaN(1:S.N);
nu_half(2:N-1)=0.5*(nu(2:S.N-1)+nu(1:S.N-2));
nu_half(1)=nu_left;
nu_half(N)=nu_right;
%Fluxes for nu equation
nu_flux=zeros(S.N,1);
nu_flux(2:S.N-1)=0.5*(u(1:S.N-2)+u(2:S.N-1));
du_dh_left=-P_L(S,nu_left)*nu_left/zeta(nu_left);
du_dh_right=-P_L(S,nu_right)*nu_right/zeta(nu_right);
nu_flux(1)=-3*du_dh_left*dh(1)/8+9*u(1)/8-u(2)/8;
nu_flux(N)=3*du_dh_right*dh(N-1)/8-9*u(N-1)/8+u(N-2)/8;
%Fluxes for u equation
u_grad=2*(u(2:S.N-2)-u(1:S.N-2))./(dh(2:S.N-1)+dh(1:S.N-2));
u_flux=NaN(S.N,1);
u_flux(2:S.N-1)=a_1*P_L(S,nu_half(2:S.N-1))+(a_2*zeta(nu_half(2:S.N-1))./nu_half(2:S.N-1)).*u_grad;
% The system of equationsY
dnu_dt = nu_flux(2:S.N)-nu_flux(1:S.N-1);
du_dt =u_flux(2:S.N)-u_flux(1:S.N-1);
% Return the derivatives
dydt = [dnu_dt; du_dt];
end
function y=P_L(S,nu)
y=S.K./nu;
end
function y=zeta(nu)
y=(2/3)*(nu).^(-2)./(nu-1);
end anonymous, function MATLAB Answers — New Questions
