Making a PINN model for adsorption column

[Rahulx084]

PINN( Carbon Capture).pdf

Please see the pdf for the partial differential equation I'm trying to solve; they're non-dimensionalized. I'm writing this code for the PDEs. The results are not coming right, the nature of the breakthrough curve is not what is expected to be. Please if someone can explain what I'm doing wrong.

import deepxde as dde
import numpy as np
import tensorflow as tf

geom=dde.geometry.Interval(0,1)
timedomain=dde.geometry.TimeDomain(0,3)
geomtime=dde.geometry.GeometryXTime(geom,timedomain)

e=0.56 #porosity
rp=0.325/10**3 #particle radius (m)
rho_p=1818 #particle density (kg/m3)
yo=0.15 #Mole fraction of CO2 in feed
L=0.2 #Bed Length (m)
uo=0.019 # Feed Velocity (m/s)
k=0.0137 # Rate constant (1/s)
Po=100000 # Feed Pressure (Pa)
vis=0.0000172 # Viscocity of gas (N/ms^2)
Co=6.03 # Feed Concentration (mol/m3)
b=1 # constant in Langmuir isotherm (m3/mol)
qsb=1.57 #Saturation concentration (mol/kg)
Pe=147 #Peclet number

b_new=bCo
m=((4rprpePo)/(150(1-e)visuoL))
n=(((1-e)rho_pqsb)/(eCo))
k_new=k*(L/uo)

def pde(inp,out):

out_c, out_q, out_p= out[:,0:1], out[:,1:2], out[:,2:] #c,q,p



out_c_t=dde.grad.jacobian(out,inp,i=0,j=1)    #dc/dt 
out_c_z=dde.grad.jacobian(out,inp,i=0,j=0)   #dc/dz
out_c_zz=dde.grad.hessian(out,inp,component=0,i=0,j=0) #d2c/dz2

out_q_t=dde.grad.jacobian(out,inp,i=1,j=1)  #dq/dt
out_q_z=dde.grad.jacobian(out,inp,i=1,j=0) #dq/dz

    
out_p_t=dde.grad.jacobian(out,inp,i=2,j=1) #dp/dt
out_p_z=dde.grad.jacobian(out,inp,i=2,j=0) #dp/dz
out_p_zz=dde.grad.hessian(out,inp,component=2,i=0,j=0) #d2p/dz2


u =((-m)*out_p_z) #velocity Darcy's law

d_cu_z= dde.grad.jacobian(out_c*u,inp,j=0) # d(cu)/dz

d_pu_z= dde.grad.jacobian(out_p*u,inp,j=0) #d(Pu)/dz



qe= ((b_new*out_c)/(1+b_new*out_c)) #Isotherm








r1= out_c_t + d_cu_z + n*out_q_t - (1/Pe)*out_c_zz

r2= out_q_t - k_new*qe + k_new*out_q

r3= out_p_t + d_pu_z + yo*n*out_q_t - (1/Pe)*out_p_zz

return [r1,r2,r3]

def boundary_l(x, on_boundary):
return on_boundary and np.isclose(x[0], 0)

def boundary_r(x, on_boundary):
return on_boundary and np.isclose(x[0],1)

bc1 = dde.icbc.DirichletBC(geomtime,lambda x :1 ,boundary_l,component=0 ) # C(0,t)=Co , 1 in function because its C'=C/Co

bc2=dde.icbc.OperatorBC(geomtime,lambda x, y, _ :(1/Pe)*dde.grad.jacobian(y, x, i=0, j=0)- y[:,0:1] + 1, boundary_l) # Danckwerts Boundary condition

bc3 = dde.icbc.DirichletBC(geomtime,lambda x :1 ,boundary_l,component=2 ) #P(0,t)=Po

bc4=dde.icbc.OperatorBC(geomtime,lambda x, y, _ : dde.grad.jacobian(y,x , i=0, j=0), boundary_r) # dc/dz(L,t)=0

ic1 = dde.icbc.IC(geomtime, lambda x: 0, lambda _, on_initial: on_initial,component=0) # C(z,0)=0

ic2 = dde.icbc.IC(geomtime, lambda x: 0, lambda _, on_initial: on_initial,component=1) # q(z,0) =0

data = dde.data.TimePDE(geomtime, pde, [bc1,bc2,bc3,bc4,ic1,ic2],
num_domain=1560, num_boundary=580, num_initial=150,num_test=100)

net = dde.nn.FNN([2] + 3*[32] + [3], "tanh", "Glorot normal")

model =dde.Model(data, net)
model.compile("adam", lr=1e-3 )
losshistory, train_state = model.train(iterations=30000)