actually I am struggling with a linear inhomogeneous differential equation of order 4 that is given by
$$ b_0\cdot\varepsilon+\sum_{i=1}^{N}b_i\frac{\partial^i \varepsilon}{\partial t^i} = a_0\cdot \sigma, \text{and } N = 1\dots 4 $$
$\sigma$ is a known constant (loading) and $a_i$, $b_i$ are functions by $t$. A - for my purpose - sufficient result can be reached when I solve the equation incrementally as IVP under the assumption that the coefficients a_i and b_i remain constant during a single time-step. This leads to the general solution
$$ \varepsilon(t) = \sum_{i=1}^{N} C_i\cdot e^{-\lambda_i t} + \frac{a_0}{b_0}\sigma $$
$\lambda_i$ are the roots of the characteristic polynomial; the constants $C_i$ are determined while solving the IVP what leads to 4x4 linear equation system. Unfortunately with highly increased time this equation system gets badly scaled and becomes unsolveable. This is just a reason of bad scaling and/or storage overflow of the variables (float, single-precision 32-bit). (In general there must exist a solution of the equation.)
Because I don't know if I manage to proper rescale the equation system, I thought about that it could be better to solve the given problem numerically by classical 4th-order RUNGE-KUTTA method. This seems also to be the better solution in terms of numbers of calculation operations.
So my questions:
- Can it be said in general that a numerical method (like RK) is a better solution in my case?
- I have to implement the method by myself in a basic-like-scripting language and therefor I would ask if there exist some open implementations that I can use or translate to my language? (Especially when I think of step-size control etc.)
- Is the "classical" RK-method also suitable for stiff IVPs? If not, do there exist some well documented literature that I can use for instance as guideline to code lets say a ROSENBROCK-WANNER-method?
- And something "general": is it possible that an ODE is stiff, even if their coefficients are constant?
Thank you in advance for your help,