# Alternatives to FiniteElement as Spatial Discretization Method for NDSolve

[...] It is possible to skip this section and continue with the discretization stage and make use of the initialized data structures ProcessEquations creates. With this it is possible to use ProcessEquations as an equation preprocessor, for example, for a new numerical discretization method.

[...]

Currently, the only discretization method available in this framework is the finite element method. Thus, by default, InitializePDEMethodData generates a FEMMethodData object.

## 1

It is possible to implement a new spatial discretization method to plug in the NDSolve framework? Does anyone tried/succeded in that? There is some reference/tutorial/skeleton implementation like the ones for time integration?

## 2

Second, related, question. Before the advent of Finite Element Framework in Mathematica 10, how Mathematica solved steady-state, boundary value problems? It is still possible to use the old method?

I know that, under NDSolve's Method option, "PDEDiscretization" only supports "MethodOfLines" and "FiniteElement". With "MethodOfLines" we can apparently enable something like a finite difference discretization with "SpatialDiscretization" -> "TensorProductGrid" but the whole "MethodOfLines" it's unsuitable for steady state, elliptic problems.

• For the second part of question, in that era, Mathematica just can not handle this type of problem, at least cannot handle it elegantly. A typical answer is this. In fact this is one of the reasons causing my self-learning of FDM and FEM. – xzczd Dec 25 '14 at 2:44
• For some more examples of how to use functionality that existed in previous versions, you may also want to look at FiniteDifferenceDerivative which I used in my answer to Finding the eigenfunctions of one and two dimensional Harmonic Oscillator and the previous post linked therein. Since much of the necessary functionality was missing, I ended up hand-rolling things like the relaxation method, e.g., here: Poisson solver using Mathematica – Jens Dec 25 '14 at 4:35
• In the past I solved some simple, 1-D, boundary value problem, with NDSolve. @user21 in his answer to my question wrote about a pre-V10 pre-FEM behavior. What's this behavior? It's still possible to use it? – unlikely Dec 25 '14 at 8:50
• If user21 doesn't appear under your question, then @ won't remind him… and I think that behavior is irrelevant to this problem. – xzczd Dec 31 '14 at 3:09
• What you can do it use NDSolve as a preprocessor, that means parse the PDE coefficients and boundary conditions. See the Finite Element programming tutorial. You'd then write your own DiscretizePDE and DiscretizeBoundaryConditions function and use LinearSolve directly. There is, however, no way to link them into NDSolve as can be done for the time integration. Would such a feature be helpful for you? If enough people request this I could think about making it available. – user21 Jan 14 '15 at 8:14