PDEToolbox - Unsuitable initial guess U0 (default: U0=0)

6 次查看（过去 30 天）
Gabriel S
Gabriel S 2020-3-29
Hello everyone.
I am trying to model a reaction of ammonia cracking ocurring in a Plug Flow (Cylindrical) Reactor with the geometry in the annex figure AQ4.jpg. For that end, I created the code below:
%1. Reactor Configuration [m]
L_aq4 = 1.5; %Reactor length
D_aq4 = 0.1; %Reactor diameter
%2. Miscellaneous Parameters
Diff_H2 = (5*10^-9)*800^1.682; %Average diffusion coefficient
lambda_gas = 6.09*10^-2; %Average thermal condutivity of the gas
%3. Creating the PDE Model
N = 5;
model = createpde(N);
% 1 = Total matter balance for the gases
% 2 = Matter balance for species NH3
% 3 = Matter balance for species N2
% 4 = Matter balance for species H2
% 5 = Heat balance
%3.2.1. Creating the PDE Model - Geometry (Cylinder)
gm = multicylinder (D_aq4/2,L_aq4);
model.Geometry = gm;
mesh_HminHmax = generateMesh(model,'Hmax',0.5,'Hmin',D_aq4);
pdegplot(model,'FaceLabels','on','FaceAlpha',0.5);
%3.2.2. Setting Boundary Conditions (Initial conditions not needed
%for steady-state problems)
C0_NH3 = 160; %[mol.m^-3]
T0 = 300; %[K]
Tw = 800; %[K]
h1 = ones(N,1);
h1d = diag(h1);
h3 = [0 0 0 0 1];
h3d = diag(h3);
r1 = [C0_NH3 C0_NH3 0 0 T0];
r3 = [0 0 0 0 Tw];
applyBoundaryCondition(model,'dirichlet','Face',1,'h',h1d,'r',r1);
applyBoundaryCondition(model,'dirichlet','Face',3,'h',h3d,'r',r3);
%3.2.3. Setting PDE Coefficients (Elliptic equation type)
c = [0 Diff_H2 Diff_H2 Diff_H2 lambda_gas]';
specifyCoefficients(model,'m',0,'d',0,'c',c,'a',0,'f',@fcoeffunction);
result = solvepde(model);
pdeplot3D(model,'ColorMapData',result.NodalSolution)
%4. Functions
function f = fcoeffunction(location,state)
nr = length(location.x);
f = zeros(5,nr);
R = 8.314; %Gas constant
ug = 1.5; %Gas velocity
cpavg_g = 3.727*10^1; %Average heat capacity of the gas
dHr_cr = 1.611*10^5; %Heat of cracking reaction
k_cr = 8.4*10^8; %Pre-exponential factor of cracking reaction
E_cr = 1.611*10^5; %Activation energy of cracking reaction
K_cr = k_cr.*exp(-E_cr./(R.*state.u(5,:))); %Arrhenius "constant" of cracking reaction
r_cr = K_cr.*state.u(2,:); %Reaction rate
f(1,:) = -ug.*state.uz(1,:) + 2.*r_cr;
f(2,:) = ug.*state.uz(2,:) + 2.*r_cr;
f(3,:) = ug.*state.uz(3,:) - r_cr;
f(4,:) = ug.*state.uz(4,:) - (3.*r_cr);
f(5,:) = (ug.*state.u(1,:).*cpavg_g).*state.uz(5,:) + dHr_cr.*r_cr;
end
The system of PDEs consist in the following equations:
1) ug*du(1)/dz = 2*r_cr
2) ug*du(2)/dz = -2*r_cr + Diff_H2*div(grad(u(2)))
3) ug*du(3)/dz = r_cr + Diff_H2*div(grad(u(3)))
4) ug*du(4)/dz = 3*r_cr + Diff_H2*div(grad(u(4)))
5) ug*cpavg_g*u(1)*du(5)/dz = -dHr_cr*r_cr + lambda*div(grad(u(5)))
For the boundary conditions:
a) For Face 1 (entrance of reactor): u(1) = u(2) = 160; u(3) = u(4) = 0; u(5) = 300;
b) For Face 3 (walls): u(5) = 800. Other variables are not bounded.
c) For Face 2 (exit): No boundary conditions.
But when I try to run the code, I get the following errors:
Warning: Matrix is singular to working precision.
> In pde.EquationModel/solveStationaryNonlinear (line 21)
In pde.PDEModel/solvepde (line 77)
In Codigo_AQ4 (line 55)
Error using pde.EquationModel/solveStationaryNonlinear (line 27)
Unsuitable initial guess U0 (default: U0=0).
By reading the documentation, the error "Unsuitable initial guess U0" means the initial guess is NaN or Inf. Is there any way to change the default U0 = 0 to prevent this? For what I understand, since I specified the boundary conditions at the entrace (Face 1), I would have specified U0, which, in this case, is not zero. Is there anything I am overlooking? Thanks in advance for the help.

请先登录，再回答此问题。

回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Structural Mechanics 的更多信息

标签

产品


版本

R2019a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by