#include <stdio.
h>
#include <stdlib.h>
// To use M_PI on windows
#define _USE_MATH_DEFINES
#include "math.h"
TwoMassesData* init_dirdyn(double ti, double tf, double dt)
{
TwoMassesData *data;
data = malloc(sizeof(TwoMassesData));
data->ti = ti;
data->tf = tf;
data->dt = dt;
data->n_steps =
data->m1 = 29;
data->m2 = 305;
data->k1 = 200;
data->k2 = 38;
data->d1 = 100;
data->d2 = 3000;
data->z1 = 0.375;
data->z2 = 0.8;
data->process = 0;
data->Fmax = 10*10^3;
data->f0 = 1;
data->t0 = 0;
data->f1 = 10;
data->t1 = 10;
TwoMassesResults *results;
results = malloc(sizeof(TwoMassesResults));
void run_dirdyn(TwoMassesData* two_masses_data) {
// A COMPLETER
void finish_dirdyn(TwoMassesData* two_masses_data) {
// A COMPLETER
void rk4(double y[], double dydx[], double x, double h, double yout[], void(*derivs)(double, double[],
double[], TwoMassesData*), TwoMassesData* s)
{
int i;
double xh, hh, h6, *dym, *dyt, *yt;
int n = 4;
dym = (double*)calloc(4, sizeof(double));
dyt = (double*)calloc(4, sizeof(double));
yt = (double*)calloc(4, sizeof(double));
hh = h*0.5;
h6 = h / 6.0;
xh = x + hh;
// First step
(*derivs)(x, yt, dydx, s);
for (i = 0; i < n; i++)
yt[i] = y[i] + hh*dydx[i];
// Second step
(*derivs)(xh, yt, dyt, s);
for (i = 0; i < n; i++)
yt[i] = y[i] + hh*dyt[i];
// Third step
(*derivs)(xh, yt, dym, s);
for (i = 0; i < n; i++)
{
yt[i] = y[i] + h*dym[i];
dym[i] += dyt[i];
}
// Fourth step
(*derivs)(x + h, yt, dyt, s);
// Accumulate increments with proper weights
for (i = 0; i < n; i++)
yout[i] = y[i] + h6*(dydx[i] + dyt[i] + 2.0*dym[i]);
// Update derivative state vector at the end of time step
(*derivs)(x, yout, dydx, s);
free(dym);
free(dyt);
free(yt);
}