diff --git a/src/ODE/Make/files b/src/ODE/Make/files
index 765fb7b207..98bbffc29c 100644
--- a/src/ODE/Make/files
+++ b/src/ODE/Make/files
@@ -2,7 +2,9 @@ ODESolvers/ODESolver/ODESolver.C
ODESolvers/ODESolver/ODESolverNew.C
ODESolvers/adaptiveSolver/adaptiveSolver.C
+ODESolvers/Euler/Euler.C
ODESolvers/EulerSI/EulerSI.C
+ODESolvers/Trapezoid/Trapezoid.C
ODESolvers/RKF45/RKF45.C
ODESolvers/RKCK45/RKCK45.C
ODESolvers/RKDP45/RKDP45.C
diff --git a/src/ODE/ODESolvers/Euler/Euler.C b/src/ODE/ODESolvers/Euler/Euler.C
new file mode 100644
index 0000000000..1a6edf576c
--- /dev/null
+++ b/src/ODE/ODESolvers/Euler/Euler.C
@@ -0,0 +1,88 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+\*---------------------------------------------------------------------------*/
+
+#include "Euler.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(Euler, 0);
+ addToRunTimeSelectionTable(ODESolver, Euler, dictionary);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::Euler::Euler(const ODESystem& ode, const dictionary& dict)
+:
+ ODESolver(ode, dict),
+ adaptiveSolver(ode, dict),
+ err_(n_)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::scalar Foam::Euler::solve
+(
+ const ODESystem& ode,
+ const scalar x0,
+ const scalarField& y0,
+ const scalarField& dydx0,
+ const scalar dx,
+ scalarField& y
+) const
+{
+ // Calculate error estimate from the change in state:
+ forAll(err_, i)
+ {
+ err_[i] = dx*dydx0[i];
+ }
+
+ // Update the state
+ forAll(y, i)
+ {
+ y[i] = y0[i] + err_[i];
+ }
+
+ return normalizeError(y0, y, err_);
+}
+
+
+void Foam::Euler::solve
+(
+ const ODESystem& odes,
+ scalar& x,
+ scalarField& y,
+ scalar& dxTry
+) const
+{
+ adaptiveSolver::solve(odes, x, y, dxTry);
+}
+
+
+// ************************************************************************* //
diff --git a/src/ODE/ODESolvers/Euler/Euler.H b/src/ODE/ODESolvers/Euler/Euler.H
new file mode 100644
index 0000000000..90f67f4b4c
--- /dev/null
+++ b/src/ODE/ODESolvers/Euler/Euler.H
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+Class
+ Foam::Euler
+
+Description
+ Euler ODE solver of order (0)1.
+
+ The method calculates the new state from:
+ \f[
+ y_{n+1} = y_n + \delta_x f
+ \f]
+ The error is estimated directly from the change in the solution,
+ i.e. the difference between the 0th and 1st order solutions:
+ \f[
+ err_{n+1} = y_{n+1} - y_n
+ \f]
+
+SourceFiles
+ Euler.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Euler_H
+#define Euler_H
+
+#include "ODESolver.H"
+#include "adaptiveSolver.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class Euler Declaration
+\*---------------------------------------------------------------------------*/
+
+class Euler
+:
+ public ODESolver,
+ public adaptiveSolver
+{
+ // Private data
+
+ mutable scalarField err_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("Euler");
+
+
+ // Constructors
+
+ //- Construct from ODE
+ Euler(const ODESystem& ode, const dictionary& dict);
+
+
+ // Member Functions
+
+ //- Solve a single step dx and return the error
+ scalar solve
+ (
+ const ODESystem& ode,
+ const scalar x0,
+ const scalarField& y0,
+ const scalarField& dydx0,
+ const scalar dx,
+ scalarField& y
+ ) const;
+
+ //- Solve the ODE system and the update the state
+ void solve
+ (
+ const ODESystem& ode,
+ scalar& x,
+ scalarField& y,
+ scalar& dxTry
+ ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/ODE/ODESolvers/Trapezoid/Trapezoid.C b/src/ODE/ODESolvers/Trapezoid/Trapezoid.C
new file mode 100644
index 0000000000..8a2c281433
--- /dev/null
+++ b/src/ODE/ODESolvers/Trapezoid/Trapezoid.C
@@ -0,0 +1,93 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+\*---------------------------------------------------------------------------*/
+
+#include "Trapezoid.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(Trapezoid, 0);
+ addToRunTimeSelectionTable(ODESolver, Trapezoid, dictionary);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::Trapezoid::Trapezoid(const ODESystem& ode, const dictionary& dict)
+:
+ ODESolver(ode, dict),
+ adaptiveSolver(ode, dict),
+ err_(n_)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::scalar Foam::Trapezoid::solve
+(
+ const ODESystem& ode,
+ const scalar x0,
+ const scalarField& y0,
+ const scalarField& dydx0,
+ const scalar dx,
+ scalarField& y
+) const
+{
+ // Predict the state using 1st-order Trapezoid method
+ forAll(y, i)
+ {
+ y[i] = y0[i] + dx*dydx0[i];
+ }
+
+ // Evaluate the system for the predicted state
+ ode.derivatives(x0 + dx, y, err_);
+
+ // Update the state as the average between the prediction and the correction
+ // and estimate the error from the difference
+ forAll(y, i)
+ {
+ y[i] = y0[i] + 0.5*dx*(dydx0[i] + err_[i]);
+ err_[i] = 0.5*dx*(err_[i] - dydx0[i]);
+ }
+
+ return normalizeError(y0, y, err_);
+}
+
+
+void Foam::Trapezoid::solve
+(
+ const ODESystem& odes,
+ scalar& x,
+ scalarField& y,
+ scalar& dxTry
+) const
+{
+ adaptiveSolver::solve(odes, x, y, dxTry);
+}
+
+
+// ************************************************************************* //
diff --git a/src/ODE/ODESolvers/Trapezoid/Trapezoid.H b/src/ODE/ODESolvers/Trapezoid/Trapezoid.H
new file mode 100644
index 0000000000..327dbe44bc
--- /dev/null
+++ b/src/ODE/ODESolvers/Trapezoid/Trapezoid.H
@@ -0,0 +1,104 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see .
+
+Class
+ Foam::Trapezoid
+
+Description
+ Trapezoidal ODE solver of order (1)2.
+
+SourceFiles
+ Trapezoid.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Trapezoid_H
+#define Trapezoid_H
+
+#include "ODESolver.H"
+#include "adaptiveSolver.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class Trapezoid Declaration
+\*---------------------------------------------------------------------------*/
+
+class Trapezoid
+:
+ public ODESolver,
+ public adaptiveSolver
+{
+ // Private data
+
+ mutable scalarField err_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("Trapezoid");
+
+
+ // Constructors
+
+ //- Construct from ODE
+ Trapezoid(const ODESystem& ode, const dictionary& dict);
+
+
+ // Member Functions
+
+ //- Solve a single step dx and return the error
+ scalar solve
+ (
+ const ODESystem& ode,
+ const scalar x0,
+ const scalarField& y0,
+ const scalarField& dydx0,
+ const scalar dx,
+ scalarField& y
+ ) const;
+
+ //- Solve the ODE system and the update the state
+ void solve
+ (
+ const ODESystem& ode,
+ scalar& x,
+ scalarField& y,
+ scalar& dxTry
+ ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //