Permalink
Branch: master
Find file Copy path
Fetching contributors…
Cannot retrieve contributors at this time
109 lines (64 sloc) 1.63 KB

Dynamics

Time discretisation

Verlet

xt::xtensor<double,2> u    = xt::zeros<double>({nnode, ndim});
xt::xtensor<double,2> v    = xt::zeros<double>({nnode, ndim});
xt::xtensor<double,2> a    = xt::zeros<double>({nnode, ndim});
xt::xtensor<double,2> v_n  = xt::zeros<double>({nnode, ndim});
xt::xtensor<double,2> a_n  = xt::zeros<double>({nnode, ndim});

xt::xtensor<double,2> fres = xt::zeros<double>({nnode, ndim});

GooseFEM::MatrixDiagonal M(...);

...

while ( ... )
{
  // history

  xt::noalias(v_n) = v;
  xt::noalias(a_n) = a;

  // new displacement

  xt::noalias(u) = u + dt * v + 0.5 * std::pow(dt,2.) * a;

  ...

  // new acceleration

  M.solve(fres, a);

  // new velocity

  xt::noalias(v) = v_n + .5 * dt * ( a_n + a );

}

Velocity Verlet

while ( ... )
{
  // variables & history

  xt::noalias(v_n) = v;
  xt::noalias(a_n) = a;

  // new displacement

  xt::noalias(u) = u + dt * v + 0.5 * std::pow(dt,2.) * a;

  ...

  // estimate new velocity

  xt::noalias(v) = v_n + dt * a_n;

  ...

  M.solve(fres, a);

  xt::noalias(v) = v_n + .5 * dt * ( a_n + a );

  ...

  // new velocity

  M.solve(fres, a);

  xt::noalias(v) = v_n + .5 * dt * ( a_n + a );

  ...

  // new acceleration

  M.solve(fres, a);

}

Forward Euler

xt::xtensor<double,2> u = xt::zeros<double>({nnode, ndim});
xt::xtensor<double,2> v = xt::zeros<double>({nnode, ndim});

...

while ( ... )
{
  xt::noalias(u) = u + dt * v;

  ...
}