I'm trying to solve a system of non-linear ODEs numerically. I'm not familiar with the methods available to do this, so I'm looking for some help to see if I'm pointing at the right direction.
The problem is the following.
Let $g=g(x)$ and $h=h(x)$ both be continuous functions of $x$, for $x\in\left[\underline{x},\infty\right)$, and let $h_{i}$ and $g_{i}$ denote their $i^{th}$ derivatives. I am looking for these functions $g$ and $h$ such that,
$0=F^{1}\left(x,h,h_{1},h_{2},h_{3},g,g_{1}g_{2}g_{3}\right)$
$0=F^{2}\left(x,h,h_{1},h_{2},h_{3},h_{7},g,g_{1}\right)$,
where $F^{1}(\cdot)$ and $F^{2}(\cdot)$ are real-valued non-linear functions of their arguments. I also have some boundary conditions for $g(\underline{x})$ and $h(\underline{x})$.
So far I was thinking about using finite differences to write down the problem, which would define a non-linear system for $h$ and $g$. Here's my main idea:
Discretize the state-space of $x$ into a grid of size $n$: $X=[\underline{x},\underline{x}+\triangle, \underline{x}+2\triangle,...,\bar{x}$].
Write down the derivatives $g_{i}$ and $h_{i}$ approximating them by finite differences. I would use a central difference scheme.
Get an expression for $F^{1}(\cdot)$ and $F^{2}(\cdot)$ in terms of $g(j)$ and $h(j)$, i.e. the functions evaluated at the $j^{th}$ entry of $X$.
Use a non-linear solver to find a solution to the system. This means finding the roots $(g(1), g(2), ... g(n))$ and $(h(1), h(2), ... h(n))$.
My main question is am I thinking this problem correctly? Should I use a particular non-linear routine for this problem? How would you include the restriction of $h(x)$ being concave? Finally, I'm using Matlab to program this, so maybe there's a command it would be worth using?
Thanks a lot!
If you can solve for the highest derivative you can use runge kutta methods. You can extend Runge kutta to solve higher derivative problems by introducing new variables. For example for a second derivative problem $y'' = F(y,y'...)$ use $v = y'$ and $v' = F(y,v,...)$.