# The heat equation solution use MethodOfLines.jl

**URL:** <https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331>\
**Category:** Modelling & Simulations\
**Tags:** question\
**Created:** [February 9, 2023, 10:29am UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331 "2023-02-09T10:29:35Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![Fred\_He](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fred_he/32/14298_2.png) [@Fred\_He](https://discourse.julialang.org/u/Fred_He)\
**Post date:** [February 9, 2023, 10:29am UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/1 "2023-02-09T10:29:35Z")

</div>

the solution looks weird in `https://docs.sciml.ai/MethodOfLines/stable/tutorials/heatss/`  
why its u(1, 1) = 0?  
 ![image](https://global.discourse-cdn.com/julialang/original/3X/d/9/d9d6f72bd9c2e29c57829fb2b2c0fa49fb1726d5.png)

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [February 9, 2023, 3:47pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/2 "2023-02-09T15:47:05Z")

</div>

Yeah, that’s odd. @xtalax

---

<div class="post-metadata">

**Author:** ![Sevi](https://avatars.discourse-cdn.com/v4/letter/s/c67d28/32.png) [@Sevi](https://discourse.julialang.org/u/Sevi)\
**Post date:** [February 9, 2023, 8:09pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/3 "2023-02-09T20:09:29Z")

</div>

I’m not sure if this is the intended behavior (still in the process of digging in the code), but it feels like the equations from `symbolic_discretize` show the wrong boundary conditions:

```julia
prob_sym = symbolic_discretize(pdesys, discretization)
for bc in prob_sym[1].eqs[end-3:end]
    @show bc
end

# bc = 0 ~ -u[1, 1]
# bc = 0 ~ -u[11, 1]
# bc = 0 ~ -u[1, 11]
# bc = 0 ~ -u[11, 11]

```

The last one should be `1 ~ u[11,11]` I guess. It also looks the same if one specifies it explicitly as

> **bcs**
>
> ```julia
> bcs = [u(0, y) ~ 0,
> u(1, y) ~ y,
> u(x, 0) ~ 0,
> u(x, 1) ~ x]
> 
> ```

Perhaps this is a starting point?

PS: The boundary conditions for the other 36 points look correct, it’s just the `u[11,11]` corner.

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [February 10, 2023, 12:37pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/4 "2023-02-10T12:37:09Z")

</div>

That would definitely be a bug. Please open an issue.

---

<div class="post-metadata">

**Author:** ![Sevi](https://avatars.discourse-cdn.com/v4/letter/s/c67d28/32.png) [@Sevi](https://discourse.julialang.org/u/Sevi)\
**Post date:** [February 10, 2023, 1:33pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/5 "2023-02-10T13:33:15Z")

</div>

Is it related to this issue?

> <https://github.com/SciML/MethodOfLines.jl/issues/148>
>
> At present corners are zeroed, better would be to interpolate/extrapolate them. …Easy in 2 dimensions, much harder to do correctly in higher dimensions.

If I understand it correctly, the corner values are not that important anyways… for (continuous) Dirichlet b.c. their value is known and for Neumann/mixed b.c. their value needs to be consistent (in some sense) with the adjacent parts of the boundary, which would require some extrapolation.

---

<div class="post-metadata">

**Author:** ![xtalax](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xtalax/32/35293_2.png) [@xtalax](https://discourse.julialang.org/u/xtalax)\
**Post date:** [February 10, 2023, 1:45pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/6 "2023-02-10T13:45:13Z")

</div>

Yes this is related. The corner equations have no effect on the interior solution, and their states do not appear in any other equation. As such they are currently held at 0, but I want to interpolate them using adjacent values so that they appear more in line with the rest of the solution.

---

<div class="post-metadata">

**Author:** ![Sevi](https://avatars.discourse-cdn.com/v4/letter/s/c67d28/32.png) [@Sevi](https://discourse.julialang.org/u/Sevi)\
**Post date:** [February 10, 2023, 2:09pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/7 "2023-02-10T14:09:00Z")

</div>

But this is only true for the “usual” stencils, right?

One could, in principle, come up with some stencils that use the corner points also for an equation in the interior. Not sure if this would be useful for anything though…

---

<div class="post-metadata">

**Author:** ![xtalax](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xtalax/32/35293_2.png) [@xtalax](https://discourse.julialang.org/u/xtalax)\
**Post date:** [February 10, 2023, 2:11pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/8 "2023-02-10T14:11:53Z")

</div>

In principle yes, and in that case a choice will need to be made about which equation governs the corner states - I don’t know of any discretization where this is the case but if you have literature I’d love to see it.

---

<div class="post-metadata">

**Author:** ![Sevi](https://avatars.discourse-cdn.com/v4/letter/s/c67d28/32.png) [@Sevi](https://discourse.julialang.org/u/Sevi)\
**Post date:** [February 10, 2023, 2:23pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/9 "2023-02-10T14:23:45Z")

</div>

I was mostly just curious about it, but apparently some people have been thinking about it at least theoretically

> <https://math.stackexchange.com/questions/2916234/how-to-obtain-the-9-point-laplacian-formula>

> **[GitHub - maroba/findiff: Python package for numerical derivatives and partial...](https://github.com/maroba/findiff)**
>
> Python package for numerical derivatives and partial differential equations in any number of dimensions. - GitHub - maroba/findiff: Python package for numerical derivatives and partial differential...

The Python package doesn’t really use this as far as I see, but they show an example of taking the derivatives in a 45º rotated basis, then the stencil looks like an “X”. But as I said, I’m not sure if there is any benefit to this and I haven’t seen it anywhere in the wild 😅

---

<div class="post-metadata">

**Author:** ![Fred\_He](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fred_he/32/14298_2.png) [@Fred\_He](https://discourse.julialang.org/u/Fred_He)\
**Post date:** [February 11, 2023, 10:46am UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/10 "2023-02-11T10:46:31Z")

</div>

BTW, it seems that `MethodOfLines.jl` only works when the coefficients are constant, for example:

```julia
dT/dt = -alpha*d^2T/dz^2

```

with `alpha` being constant. It will not work if the coefficients are functions of the dependent variables, for example:

```julia
dT/dt = -alpha(T)*d^2T/dz^2

```

with `alpha` being dependent on `T`

---

<div class="post-metadata">

**Author:** ![xtalax](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xtalax/32/35293_2.png) [@xtalax](https://discourse.julialang.org/u/xtalax)\
**Post date:** [February 11, 2023, 3:27pm UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/11 "2023-02-11T15:27:30Z")

</div>

> BTW, it seems that `MethodOfLines.jl` only works when the coefficients are constant

No, this should work. Did you `@register_symbolic alpha(T)`?

---

<div class="post-metadata">

**Author:** ![Fred\_He](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fred_he/32/14298_2.png) [@Fred\_He](https://discourse.julialang.org/u/Fred_He)\
**Post date:** [February 13, 2023, 1:13am UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/12 "2023-02-13T01:13:13Z")

</div>

Thanks!  
Then should I define for example  
`alpha(T) = a*T+b`  
first, then  
`@register_symbolic alpha(T)`  
, and then define the equation?  
Also, which definition is correct (or recommended) ?  
`eq = Dt(T(x,t)) ~ alpha(T(x,t))*Dxx(T(x,t))`  
or  
`eq = Dt(T) ~ alpha(T)*Dxx(T)` ?

---

<div class="post-metadata">

**Author:** ![xtalax](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/xtalax/32/35293_2.png) [@xtalax](https://discourse.julialang.org/u/xtalax)\
**Post date:** [February 13, 2023, 11:31am UTC](https://discourse.julialang.org/t/the-heat-equation-solution-use-methodoflines-jl/94331/13 "2023-02-13T11:31:56Z")

</div>

The first method is generally better, you want to define your `T` like `@variables T(..)` to ease defining your bcs and ics, then `T` on its own is just the operation