[soes]

Hello
I have attached a description of the project I am working on.

I don't get my expected results, I doubt that my way of implementation of the bottom boundary condition might not be correct.

[FMDenaro]

we already discussed your problem in other posts and I still do not understand how the number of BC.s you used can close your fourth-order PDE...
Then, what about a stability analysis of your discretization?

[soes]

Yes, I mean we have two initial conditions for time and two boundary conditions at top and bottom, and periodic boundary condition in the xi direction, that's how it closes the equation.

I did von neuman analysis and it showed that for the time steps and grid size that I am using it is stable, no worries for that, even I tried crank nicolson method it gave the same results.

[FMDenaro]

Accepting "per fides" that your problem is mathematically well posed, the fact that by decreasing the steps dx and dz the solution becomes wrong is possible if:

1) the dt is too large and produces a numerical instabilty. I have not checked carefully your discrete equation, since your multi-step method produces a 2x2 amplification matrix are you sure that both eigenvalues ensure a conditional stability? However, that seems not your case as you wrote that the solution is stable and dissipate the wave.

2) the discretized problem (equation + bc.s) is not consistent, that means your local truncation error does not vanish for vanishing mesh sizes.

Could you plot your numerical solutions and the analytical one at a certain time-step?

[soes]

Here I have attached the time evolution of the solution and the exact analytical expression, as you see they perfectly mach, and the other plot shows the convergence as you can see the difference goes to almost zero 1e-6, although I should have plotted it logarithmically to show it better, but it well converges and follows that analytical solution. This is for the case without bottom corrugation.

[FMDenaro]

1) from the plot I cannot see the shift of the wave from the position at t=0, could you superimpose the initial condition?

2) the plot is just for a value of w in one point, could you plot the global solution w(x,z)?

[FMDenaro]

However, I searched in the literature but I did not find any PDE resembling that in your problem...

[mprinkey]

I agree with Prof Denaro. The splitting that we discussed before assumes that you can solve a forced wave equation for phi that is defined as Laplacian(w), extract w from a Poisson equation, and then compute the next timestep driving function for the wave equation from that solution. And, this may be a favorable splitting for the boundary conditions in the initial test problem, but not for more general BCs.

There are two related and unresolved issues...characteristic/information flow directions/domains of dependence...elliptic, parabolic, etc...and well-posedness. I have no experience with equations of this sort...and your reticence to tell us about the physical system it purports to represent doesn't offer me or the others here any impetus to develop it. Quite frankly, we can write down PDEs at random that probably have very interesting and convoluted behaviors, but we have better things to think about.

Having said all of that, let me just say that you can likely get to the bottom of at least the information flow issue by building a fully implicit time discretization for w. That would entail using an implicit time integration scheme (backward Euler, Adams Moulton, or implicit RK) and apply it to the first-order split system that we outlined in previous threads. The system may be block implicit in w and dw/dt. But, this will give you a finite difference formulation that does NOT assume any directional information flow. The next timestep information will depend completely on old time and new time information. This does assume that the system is an evolution-type equation...the domain of dependence for new time lies in the old time domain and includes only THE NEXT time information...not all future time as well. That is probably a safe assumption. And this should give you results that don't rely on SPATIAL direction assumptions. If the linear system that you build is not singular or indeterminate, that (at least) gives you an indication that the system is possibly well-posed.

I will tell you that using implicit timestepping will likely introduce diffusion into the system and your waves will be damped over time. There may be energy conserving time integration schemes (especially implicit RK) that may help with this.

[soes]

https://drive.google.com/file/d/0B4b83ymv5MKZc3ZZWm5Qb01HY28/view?usp=sharing

@FMDenaro I have simulated the one with normal boundary conditions, and this is the video of it

[soes]

Originally Posted by mprinkey

I agree with Prof Denaro. The splitting that we discussed before assumes that you can solve a forced wave equation for phi that is defined as Laplacian(w), extract w from a Poisson equation, and then compute the next timestep driving function for the wave equation from that solution. And, this may be a favorable splitting for the boundary conditions in the initial test problem, but not for more general BCs.

There are two related and unresolved issues...characteristic/information flow directions/domains of dependence...elliptic, parabolic, etc...and well-posedness. I have no experience with equations of this sort...and your reticence to tell us about the physical system it purports to represent doesn't offer me or the others here any impetus to develop it. Quite frankly, we can write down PDEs at random that probably have very interesting and convoluted behaviors, but we have better things to think about.

Having said all of that, let me just say that you can likely get to the bottom of at least the information flow issue by building a fully implicit time discretization for w. That would entail using an implicit time integration scheme (backward Euler, Adams Moulton, or implicit RK) and apply it to the first-order split system that we outlined in previous threads. The system may be block implicit in w and dw/dt. But, this will give you a finite difference formulation that does NOT assume any directional information flow. The next timestep information will depend completely on old time and new time information. This does assume that the system is an evolution-type equation...the domain of dependence for new time lies in the old time domain and includes only THE NEXT time information...not all future time as well. That is probably a safe assumption. And this should give you results that don't rely on SPATIAL direction assumptions. If the linear system that you build is not singular or indeterminate, that (at least) gives you an indication that the system is possibly well-posed.

I will tell you that using implicit timestepping will likely introduce diffusion into the system and your waves will be damped over time. There may be energy conserving time integration schemes (especially implicit RK) that may help with this.

@mprinkey
Thanks for the reply, if you have a look at my uploaded file on this post you see that I did not do any splitting and I solved the whole equation. Worth mentioning even with splitting, no matter in what way, the results are the same, upon similar resolutions.
I am not reticent in giving information I have mentioned all I know, all I was assigned by my course supervisor was to solve the PDE, I do not know more that what I have written.
Besides, I as well tried crank nicolson scheme, but still I had the limitation of spatial resolution, and it is weird, when I increase the spatial resolution, say more than 11 points per wave length the wave does not move at all, and gets damped away.

[FMDenaro]

Originally Posted by Soheil.esmaeilzadeh

https://drive.google.com/file/d/0B4b83ymv5MKZc3ZZWm5Qb01HY28/view?usp=sharing

@FMDenaro I have simulated the one with normal boundary conditions, and this is the video of it

Yes, I see now, it appears nothing else that like rigid moving wave along the z direction ... it is like solving the equation dw/dt + c* dw/dz=0 at several positions along x with sinusoidal initial condition.

That increases my doubts ... can a second order in time PDE produce a single wave solution? should not be unstable (or presenting wiggles) the solution with central discretization for the case of a wave propagation?
As I wrote previously, despite my searching, I have not found nothing in literature like your PDE.

[soes]

Originally Posted by FMDenaro

Yes, I see now, it appears nothing else that like rigid moving wave along the z direction ... it is like solving the equation dw/dt + c* dw/dz=0 at several positions along x with sinusoidal initial condition.

That increases my doubts ... can a second order in time PDE produce a single wave solution? should not be unstable (or presenting wiggles) the solution with central discretization for the case of a wave propagation?
As I wrote previously, despite my searching, I have not found nothing in literature like your PDE.

you mean a central difference might not be a good way of discretizing a wave like equation like this? since it cannot capture instability?

[FMDenaro]

Originally Posted by Soheil.esmaeilzadeh

you mean a central difference might not be a good way of discretizing a wave like equation like this? since it cannot capture instability?

if you work with hyperbolic equations, the solution has a domain of dependence defined by the characteristic curves, the central discretization for a first order PDE is clearly a bad choice as it violates the domain of dependence. Conversely, if you have high order PDE you can have solution with two waves coming from different regions and the central discretization could get reasonable result.

To tell the true, I suppose that could be the particular choice of your initial condition to have the form of an eigenfunction of the operators in your PDE...

[soes]

Originally Posted by FMDenaro

if you work with hyperbolic equations, the solution has a domain of dependence defined by the characteristic curves, the central discretization for a first order PDE is clearly a bad choice as it violates the domain of dependence. Conversely, if you have high order PDE you can have solution with two waves coming from different regions and the central discretization could get reasonable result.

To tell the true, I suppose that could be the particular choice of your initial condition to have the form of an eigenfunction of the operators in your PDE...

Then I will as well try with a forward and backward discretization (depending on closeness to each side) as well to see whether it fixes the problem that happens with increase in resolution.

[FMDenaro]

Originally Posted by Soheil.esmaeilzadeh

Then I will as well try with a forward and backward discretization (depending on closeness to each side) as well to see whether it fixes the problem that happens with increase in resolution.

no, this is not the issue...in case of first order hyperbolic equations, forward/backward stencils are automaticaly determined by the sign of the eigenvalues of the PDE, they are not determined by the limit of the domain.

[soes]

Originally Posted by FMDenaro

no, this is not the issue...in case of first order hyperbolic equations, forward/backward stencils are automaticaly determined by the sign of the eigenvalues of the PDE, they are not determined by the limit of the domain.

Then it comes again the issue, that my equation is not any of the conventional PDE types. Then I should probably switch to finite element method for instance

[FMDenaro]

Originally Posted by Soheil.esmaeilzadeh

Then it comes again the issue, that my equation is not any of the conventional PDE types. Then I should probably switch to finite element method for instance

no, finite element is just a different way to discretize the PDE, does not change the original nature of your PDE... what is more, the FE require to produce a suitable functional form of your PDE with proper shape function (in space). And the time integration is often nothing else that a FD discretization.

[soes]

Originally Posted by FMDenaro

no, finite element is just a different way to discretize the PDE, does not change the original nature of your PDE... what is more, the FE require to produce a suitable functional form of your PDE with proper shape function (in space). And the time integration is often nothing else that a FD discretization.

But a colleague has successfully with high resolution with finite element method has solved this PDE with no such issues that I have, and he has not examined the stability or has not found out the PDE type and etc, just has discretized and solved it.

[mprinkey]

The formulation in the pdf at the top of this thread is NOT implicit in time even though you are solving a linear system for the new time values of w. You used central differencing of the d^2/dt^2 operator, so the new time value is known explicitly by solving for phi_n+1 = laplacian(w_n+1). The fact that you have implicitly discretized the spatial operator and built a linear system does not make this an implicit time stepping scheme. So, this assumes some information flow in the system, just like the earlier splitting approach did.

My recommendations from before stand. Adams Moulton, backward Euler, or Implicit RK for the time advancement.

[soes]

Originally Posted by mprinkey

The formulation in the pdf at the top of this thread is NOT implicit in time even though you are solving a linear system for the new time values of w. You used central differencing of the d^2/dt^2 operator, so the new time value is known explicitly by solving for phi_n+1 = laplacian(w_n+1). The fact that you have implicitly discretized the spatial operator and built a linear system does not make this an implicit time stepping scheme. So, this assumes some information flow in the system, just like the earlier splitting approach did.

My recommendations from before stand. Adams Moulton, backward Euler, or Implicit RK for the time advancement.

Yes in the pdf there is an explicit scheme, but I have also tried crank nicolson, but it had a similar behavior, limitations and problems described above.