2
$\begingroup$

As per suggestions of the people from MathOverflow, I'm reposting my question here:

I'm currently trying to solve a linear system $Ax = B$, where the matrix $A$ is ill conditioned (i.e. nearly singular), with a condition number of $~10^7$. The aforementioned linear system arises from a finite difference discretization.

The mathematical model for my problem is a PDE with derivatives of $x$ and $t$. Therefore, I'm solving a linear system of a discretized mesh of points with interval $\Delta x$ for each time step $\Delta t$. I'm already using centered differences for time, but because the problem has mixed derivatives (i.e. dependent of both $x$ and $t$), the problem still falls into a linear system.

My question is: how can I find a solution to this linear system? The most common solution I saw is doing preconditioning on matrix $A$, so it gets better conditioned. But, due to my engineering background, I hardly understand what must be done to precondition matrix $A$. I saw many different methods (like Jacobi and ILU factorization), but I don't know how to apply them.

On a sidenote, I'm doing this on MATLAB, so if anyone know any in-built function that can help me, it would be very appreciated. I tried gmres but it didn't work (the solution still "exploded" after a few steps).

EDIT: Here is the equation as asked

$$ EIr\left[\frac{\partial^4\theta}{\partial x^4}-6\left(\frac{\partial\theta}{\partial x}\right)^2\frac{\partial^2\theta}{\partial x^2}\right]-EAr\left[\frac{\partial u_x}{\partial x}\frac{\partial^2\theta}{\partial x^2}+\frac{\partial^2 u_x}{\partial x^2}\frac{\partial\theta}{\partial x}+1.5r^2\frac{\partial^2\theta}{\partial x^2}\left(\frac{\partial\theta}{\partial x}\right)^2\right]-I_pr\omega\left[2\frac{\partial^2\theta}{\partial x \partial t}\frac{\partial\theta}{\partial x}+\frac{\partial\theta}{\partial t}\frac{\partial^2\theta}{\partial x^2}\right]-m_pr\frac{\partial u_x}{\partial x}\frac{\partial^2\theta}{\partial t^2}+m_pr\frac{\partial^2 u_x}{\partial t^2}\frac{\partial\theta}{\partial x}-(x+u_x)m_pr\frac{\partial^3\theta}{\partial x \partial^2 t}+(x+u_x)m_pr\left(\frac{\partial\theta}{\partial t}\right)^2\frac{\partial\theta}{\partial x}+m_pg\sin\theta = 0 $$

$\endgroup$
  • 1
    $\begingroup$ Are you sure that your discretization is producing a stable iteration? If it isn't, then it won't matter how accurately you solve Ax=b, you'll still be in trouble. If you back up and discuss the discretization of the problem in more detail, we might be able to answer the question at that level. $\endgroup$ – Brian Borchers Mar 5 '16 at 19:11
  • $\begingroup$ Are you perhaps refering to the selection of steps deltax and deltat? If so, yes, they respect the CFL condition for the wave equation. I used centered Euler approximations for both space and time (i.e. only dependant of the two adjacent points). $\endgroup$ – Martios Mar 5 '16 at 20:53
  • $\begingroup$ How many spatial unknowns does the system have? Depending on the problem size, you may be better off with sparse direct methods or, as you mention, an iterative method like GMRES. $\endgroup$ – Daniel Shapero Mar 5 '16 at 22:01
  • $\begingroup$ I created a mesh of 101 points. I tried using GMRES, but the solution still blew up. $\endgroup$ – Martios Mar 5 '16 at 23:02
  • 1
    $\begingroup$ Martios, you can type the equation in the question using MathJax. That way every person can see it. $\endgroup$ – nicoguaro Mar 6 '16 at 18:43
1
$\begingroup$

The two simplest solutions to solving ill posed problems are 1- truncated SVD, and 2 - Tikhonov regularization, which are actually related to one another but implemented in different ways. Plenty of details are available online about both of these methods; Wikipedia should be sufficient to get you going. I would start with Tikhonov regularization, and then move to #1 if you want a fully automatic method.

$\endgroup$
  • $\begingroup$ I see, I will take a look, thanks. Does MATLAB has any in-built function for this? $\endgroup$ – Martios Mar 6 '16 at 16:53

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.