#1-D Heat equation/diffusion equation with normalization constraint

1 messages · Page 1 of 1 (latest)

onyx sable
#

$$
\frac{du(t,x)}{dt} = D \frac{d^2 u(t,x)}{dx^2} \\
$$
$$
u(t,0) = z(t,0) \\
$$
$$
u(0,x) = 0.
$$

$$
z(t,x) = \int_0^{x} u(t,u) du \\
$$
$$
\frac{dz(t,x)}{dx} = u(t,x) \\
$$
$$
z(t,+\infty) = \phi_0 \\
$$
$$
z(0,x) = 0
$$

boreal quarryBOT
#

LonelyTree

lost lark
#

what's the finalloss?

onyx sable
# lost lark what's the finalloss?

My boundary condition was zero everywhere, that was the problem. My issue now is that I want to impose non-negativity on the concentration, how can I do this?

My updated code and equations:

function dirac(t,x)
    if isapprox(t,zero(t)) && isapprox(x,zero(x))
        return 1.
    end
    return 0.
end
@register_symbolic dirac(t::Real,x::Real)

function problem(D,ϕ_0,t_span=(0.0,120*60.0),x_span=(0.0,100.0);dt=10*60.,dx=5.,lifespan=10)
    @parameters x
    @variables t u(..) z(..)
    Dxx = Differential(x)^2
    Dx = Differential(x)
    Dt = Differential(t)

    eq  = [
        Dt(u(t,x)) ~ D*Dxx(u(t,x)),
        Dx(z(t,x)) ~ u(t,x)
    ]

    # Initial and boundary conditions
    bcs = [
        u(0,x) ~ ϕ_0*dirac(0,x),
        z(t,x_span[2]) ~ ϕ_0,
        z(t,x_span[1]) ~ 0.,
    ]

    # Space and time domains
    domains = [t ∈ t_span,
            x ∈ x_span]

    @named pde_system = PDESystem(eq,bcs,domains,[t,x],[u(t,x),z(t,x)])
    # Method of lines discretization 
    order = 2
    discretization = MOLFiniteDifference([x => dx], t)
    prob = discretize(pde_system,discretization)
    
    # Convert the PDE problem into an ODE problem
    sol = solve(prob, Rosenbrock23(), saveat=dt)

    sol, u, x, z, t
end
#

$$
\frac{du(t,x)}{dt} = D \frac{d^2 u(t,x)}{dx^2} \\
$$
$$
u(0,x) = \phi_0 \delta(0,x)
$$

$$
z(t,x) = \int_0^{x} u(t,u) du \\
$$
$$
\frac{dz(t,x)}{dx} = u(t,x) \\
$$
$$
z(t,+\infty) = \phi_0 \\
$$
$$
z(t,0) = 0
$$

boreal quarryBOT
#

LonelyTree

onyx sable
#

This is the behaviour of the normalization constraint z(t,x)

#

sol, u, x, z, t = problem(0.1,1.,(0,1),(0,1);dt=0.1,dx=0.1);