Boundary conditions for variable length bar

[hunt_mat]

TL;DR Summary

How can you supply the boundary condition for the conservation of mass to a bar of time-dependent length?

Suppose I'm looking at a bar of length $L(t)$ in 1D and I have the conservation of mass:
$$\frac{\partial\rho}{\partial t}+\frac{\partial}{\partial x}(\rho u)=0$$
In order to make things easier, I make the change of variable $x'=x/L(t)$ so that in this frame of reference, the length remains constant, and it will keep the numerics easier. The equation is then transformed into:
$$L^{2}(t)\frac{\partial\rho}{\partial t}-L'(t)\frac{\partial\rho}{\partial x'}+\frac{\partial}{\partial x'}(\rho u)=0$$
Now when I do a method of lines numerical method I end up with the following:
$$\frac{d\rho_{i}}{dt}=\frac{L'}{2hL^{2}}(\rho_{i+1}-\rho_{i-1})-\frac{L}{2h}(\rho_{i+1}u_{i+1}-\rho_{i-1}u_{i-1})$$
I thought about using the boundary condition:
$$\frac{\partial}{\partial x}(\rho u)=0$$
but I don't know how to deal with $\partial\rho/\partial x'$ on the boundary. Does anyone know?

 

[Paul Colby]

Isn't $\rho = 0$ at the boundary? That's where your bar ends.

 

[hunt_mat]

Wouldn't that simply be the value of $\rho$ at $x'=0^{-},1^{+}$?

 

[Paul Colby]

I think so. Wouldn't your difference equation then be solvable given an initial mass and velocity distribution at $t=0$?

I'm assuming there is no bar outside of the interval $[0,1]$

 

[hunt_mat]

I'm a little confused as to what you're saying, to put everything into ode45, you need to include the terms at $i=1$ and $i=N$ which means that you have to use ghost cells $\rho_{0}$ and $\rho_{N+1}$ which you normally find from the boundary conditions, I don't kbnow how to deal with $\rho$ by itself.

 

[Paul Colby]

Let's pull back a bit. Both $\rho$ and $u$ depend on both $t$ and $x$. What you start with is an equation of conservation. Given $u(t,x)$ for all time one could in principle compute $\rho(t,x)$, which is where the mass density moves to by solving the differential equation. I don't see any equations that determine $u$?

 

[Paul Colby]

Let me spell out my issue in a little more detail. Make both $t$ and $x$ discrete. So the continuity equation is,
$$
\rho(n+1,m) - \rho(n,m) = (\rho(n,m+1)-\rho(n,m))u(n,m) + \rho(n,m)(u(n,m+1)-u(n,m))
$$
Everything in this equation is assumed given at $t=n$. So we can compute $\rho(n+1,m)$. Sweet. Now what determines $u(n+1,m)$? That seems to be missing independent of any question regarding boundary conditions.

 

[pasmith]

TL;DR Summary: How can you supply the boundary condition for the conservation of mass to a bar of time-dependent length?

Suppose I'm looking at a bar of length $L(t)$ in 1D and I have the conservation of mass:
$$\frac{\partial\rho}{\partial t}+\frac{\partial}{\partial x}(\rho u)=0$$
In order to make things easier, I make the change of variable $x'=x/L(t)$ so that in this frame of reference, the length remains constant, and it will keep the numerics easier. The equation is then transformed into:
$$L^{2}(t)\frac{\partial\rho}{t}-L'(t)\frac{\partial\rho}{\partial x'}+\frac{\partial}{\partial x'}(\rho u)=0$$

I think this should be
$$\frac{\partial \rho}{\partial t} - \frac{x'L'(t)}{L} \frac{\partial \rho}{\partial x'} + \frac{1}{L}\frac{\partial}{\partial x'} (\rho u) = 0$$

Now when I do a method of lines numerical method I end up with the following:
$$\frac{d\rho_{i}}{dt}=\frac{L'}{2hL^{2}}(\rho_{i+1}-\rho_{i-1})-\frac{L}{2h}(\rho_{i+1}u_{i+1}-\rho_{i-1}u_{i-1})$$
I thought about using the boundary condition:
$$\frac{\partial}{\partial x}(\rho u)=0$$
but I don't know how to deal with $\partial\rho/\partial x'$ on the boundary. Does anyone know?

On the boundary you are trying to enforce this condition: $$\frac{\partial \rho}{\partial t} - \frac{x'L'}{L}\frac{\partial \rho}{\partial x'} = 0.$$ At x = x' = 0 this is not a problem: it reduces to $\frac{\partial \rho}{\partial t} = 0$. At $x' = 1$ you can use the integrated conservation equation in the form $$\frac{d}{dt}\int_0^{L(t)} \rho\,dx = 0$$ to obtain $$L'(t) \rho(L(t)) + \int_0^L \frac{\partial \rho}{\partial t}\,dx = 0.$$ The integral vanishes because it is equal to $\left[- \rho u \right]_0^L = 0$. That leaves $$L'(t)\rho(L(t)) = 0$$ so either the bar is of constant length or $\rho(L) = 0$.

Note that you need somehow to determine $L'(t)$.

Let's pull back a bit. Both $\rho$ and $u$ depend on both $t$ and $x$. What you start with is an equation of conservation. Given $u(t,x)$ for all time one could in principle compute $\rho(t,x)$, which is where the mass density moves to by solving the differential equation. I don't see any equations that determine $u$?

Usually in these problems $u = u(\rho)$ is given. Or you can use conservation of momentum.

 

[hunt_mat]

I think this should be
$$\frac{\partial \rho}{\partial t} - \frac{x'L'(t)}{L} \frac{\partial \rho}{\partial x'} + \frac{1}{L}\frac{\partial}{\partial x'} (\rho u) = 0$$
On the boundary you are trying to enforce this condition: $$\frac{\partial \rho}{\partial t} - \frac{x'L'}{L}\frac{\partial \rho}{\partial x'} = 0.$$ At x = x' = 0 this is not a problem: it reduces to $\frac{\partial \rho}{\partial t} = 0$. At $x' = 1$ you can use the integrated conservation equation in the form $$\frac{d}{dt}\int_0^{L(t)} \rho\,dx = 0$$ to obtain $$L'(t) \rho(L(t)) + \int_0^L \frac{\partial \rho}{\partial t}\,dx = 0.$$ The integral vanishes because it is equal to $\left[- \rho u \right]_0^L = 0$. That leaves $$L'(t)\rho(L(t)) = 0$$ so either the bar is of constant length or $\rho(L) = 0$.

Note that you need somehow to determine $L'(t)$.
Usually in these problems $u = u(\rho)$ is given. Or you can use conservation of momentum.

You're correct about the error in transformation, and I agree with your analysis at $n'=0$ but I am a little confused about your analysis, I did a similar analysis for the global conservation of mass. Like you I obtain:
$$\frac{d}{dt}\int_{0}^{L(t)}\rho dx=0\Rightarrow \rho(t,L(t))L'(t)+\int_{0}^{L(t)}\frac{\partial \rho}{\partial t}dx$$
I then use the equation to obtain:
$$\rho(t,L(t))L'(t)+\int_{0}^{L(t)}\frac{\partial \rho}{\partial t}dx=\rho(t,L(t))L'(t)-\int_{0}^{L(t)}\frac{\partial}{\partial x}(\rho u)dx=\rho(t,L(t))L'(t)-(\rho(t,L(t))u(t,L(t))-\rho(t,0)u(t,0))=\rho(t,L(t))L'(t)-\rho(t,L(t))u(t,L(t))=0$$
Which yields the condition:
$$L'(t)=u(t,L(t))$$
This implicitly assumed that $\rho(t,L(t))\neq 0$ though

 

[hunt_mat]

Let's pull back a bit. Both $\rho$ and $u$ depend on both $t$ and $x$. What you start with is an equation of conservation. Given $u(t,x)$ for all time one could in principle compute $\rho(t,x)$, which is where the mass density moves to by solving the differential equation. I don't see any equations that determine $u$?

Sorry, I didn't get back to this, there are three equations for $\rho$,$u$ and $T$ but the problem appears in the conservation of mass equation, so I thought that this would be the best equation to demonstrate the problem.

 

[Paul Colby]

The equation of continuity that appears in this thread is equally valid for a solid, liquid or gaseous rod. Is it reasonable to expect the boundary conditions to be independent of this?

 

[hunt_mat]

Good question. On the free end, I would expect we can apply a stress-free condition as it's allowed to move. I'm not sure how this transfers to the density though.

 

[Paul Colby]

Well, I suggest the answer is your boundary conditions depend on the system being modeled. The question as you’ve posed here is incomplete and cannot be answered. The length of the bar changes depending on the velocity at the boundary. This is usually determined by the dynamics, which is missing.

 

[hunt_mat]

I don't think it does to be honest, the system I'm modelling is essentially the thermal poroelestic system if that helps. I'm treating it as a continuum. As a base level, I require that no mass leaves or enters the system, so the mass flux should be zero at both ends should be zero, and that's modelled by the boundary condition I gave(I think). There is a stress-free condition on the free end, so it's allowed to move.

Do you want to see the rest of the system?

 

[Paul Colby]

Not really, it’s your problem. Just in general, if I have a liquid, the manner of constraint of that liquid isn’t determined by it being a liquid or by it being conserved. The constraints, even in the case of being unconstrained, is determined by the problem. Given that your ends are unconstrained, the stress free condition says a solid will sit there oscillating while a gas will expand into a vacuum. If these are the result of your bc, you’re done.

 

[hunt_mat]

One of the ends is constrained and the other is free to move. I have a thermal poroelastic medium in which the porosity decreases as you heat up the bar. I don't need to model anything outside the bar.

 

[Paul Colby]

I don't need to model anything outside the bar.

I agree because you’ve just specified one end is fixed and one end is moving but unconstrained. The stress vanishes at this end. What does the stress equation say about $u$ and $\rho$ at this point?

 

[Frabjous]

@FEAnalyst might have some insights.