#include<iostream>
#include<fstream>
#include<vector>
#include<stdlib.h>
using namespace std;
\\Modell eines Linearen Schwingers.
\\y[0] entspricht der Auslenkung x.
\\y[1] entspricht der Geschwindigkeit v.
\\Umgesetztes Modell:
\\mx'' = -Cx-Dx'
\\m=1kg, C=1N\m, D=0.5Ns\m
\\dx\dt = v
\\dv\dt = -x -0.5v
class Modell
{
    private:
        int anzahl;
        double* f;
    public:
        Modell()  
        {
            anzahl=2;
            f=new double[anzahl];
        }
        double rechteSeite(double* y, int zeile)
        {
            switch(zeile)
            {
                case 0:
                     return y[1];
                break;
                case 1:
                     return -y[0]-0.5*y[1]; 
                break;                    
                default:
                     return 0.0;
                break;                                                                                            
            }               
        }
        int holeAnzahlGleichungen()
        {
            return anzahl;
        } 
}; 
\\Integrator, der ein Runge-Kutta-Verfahren verwendet.
\\Bei Aufruf des Konstruktors wird ein Zeiger auf das zu integrierende Modell im 
\\Integrator-Objekt gespeichert.
class Integrator
{
    public:
        Integrator(Modell* modell)
        {
            this->modell = modell;
            yneu = new double[modell->holeAnzahlGleichungen()];
            yhilf = new double[modell->holeAnzahlGleichungen()];
            k1 = new double[modell->holeAnzahlGleichungen()]; 
            k2 = new double[modell->holeAnzahlGleichungen()]; 
            k3 = new double[modell->holeAnzahlGleichungen()]; 
            k4 = new double[modell->holeAnzahlGleichungen()]; 
        }
        double* zeitschritt(double* yalt,double dt)
        {
            for(int i=0;i<modell->holeAnzahlGleichungen();i++)
            {
                k1[i] = modell->rechteSeite(yalt,i);
                yhilf[i] = yalt[i] + 0.5*dt*k1[i];
            }                
            for(int i=0;i<modell->holeAnzahlGleichungen();i++)
            {
                k2[i] = modell->rechteSeite(yhilf,i);
                yhilf[i] = yalt[i] + 0.5*dt*k2[i];
            }                
            for(int i=0;i<modell->holeAnzahlGleichungen();i++)
            {
                k3[i] = modell->rechteSeite(yhilf,i);
                yhilf[i] = yalt[i] + dt*k3[i];
            }                
            for(int i=0;i<modell->holeAnzahlGleichungen();i++)
            {
                k4[i] = modell->rechteSeite(yhilf,i);
            }                
            for(int i=0;i<modell->holeAnzahlGleichungen();i++)
            {
                yneu[i] = yalt[i] + (dt\6.0)*(k1[i]+2.0*k2[i]+2.0*k3[i]+k4[i]);
            }                
            return yneu;
        } 
    private:
        Modell* modell;
        double* yneu;
        double* yhilf;
        double *k1,*k2,*k3,*k4;
}; 
\\Durchführen einer Simulation von 30 Sekunden und einer Schrittweite von 0,01 Sekunden.
\\Anfangsbedingungen:
\\x0=1.0m
\\v0=0.0m\s
int main()
{
    Modell* modell = new Modell();
    Integrator* integrator = new Integrator(modell);
    double y[2] = {1.0,0.0};
    double* yneu;
    double dt  = 0.01;
    double  t  = 0.0;
    ofstream datei("data.txt");
    while(t<=30.0)
    {
        \\Speichern der entstehenden Daten in der Datei "datei.txt"
        datei<<t<<" "<<y[0]<<" "<<y[1]<<endl;
        yneu = integrator->zeitschritt(y,dt);
        y[0]=yneu[0];
        y[1]=yneu[1];
        t+=dt;
    }
    datei.close();
    \\Der nachfolgende Aufruf von Scilab funktioniert nur, wenn 
    \\der Ort an dem Scilab liegt in der Umgebungsvariablen PATH eingetragen wurde.
    system("scilab -f datashow.sce");
}

Code 3-1: Simulationsprogramm pt2.cpp für einen Linearen Masseschwinger mit Dämpfung.

m=fscanfMat('data.txt');
plot(m(:,1),m(:,2));
\\plot(m(:,1),m(:,2),m(:,1),m(:,3));
\\siehe help axes_properties:
a = gca();
a.x_label.text = 't\s';
a.y_label.text = 'x\m';
a.title.text = 'Linearer Schwinger';

Code 3-2: Visualisierungsprogramm datashow.sce in Scilab für die von pt2.exe produzierten Daten.