bubble-dynamics

Spherical bubble dynamics simulator
git clone https://git.0xfab.ch/bubble-dynamics.git
Log | Files | Refs | README | LICENSE

RayleighPlesset.h (2315B)


      1 /* File:   RayleighPlesset.h */
      2 /* Date:   Tue Jan 12 23:00:47 2016 */
      3 /* Author: Fabian Wermelinger */
      4 /* Tag:    Rayleigh-Plesset dynamics */
      5 /* Copyright 2016 ETH Zurich. All Rights Reserved. */
      6 #ifndef RAYLEIGHPLESSET_H_DLCI4XQ3
      7 #define RAYLEIGHPLESSET_H_DLCI4XQ3
      8 
      9 #include "BubbleBase.h"
     10 
     11 #include <ODETB/GnuplotDump.h>
     12 #include <ODETB/TimeStepper/KernelBase.h>
     13 
     14 #include <cassert>
     15 #include <cmath>
     16 
     17 template <typename Tinput, typename Trhs=Tinput>
     18 class RayleighPlesset : public KernelBase<Tinput,Trhs>
     19 {
     20     const size_t _N;
     21     GnuplotDump* m_dumper;
     22 
     23 public:
     24     RayleighPlesset(const size_t N=0, const int format=GnuplotDump::BIN) : _N(N) { m_dumper = new GnuplotDump("out", format); }
     25     virtual ~RayleighPlesset() { delete m_dumper;}
     26 
     27     static void IC(Tinput& U, BubbleData const& data)
     28     {
     29         assert(Tinput::DataType::SIZE == 2);
     30         assert(data.Nbubbles == 1);
     31         typename Tinput::DataType& IC = U[0];
     32         IC[0] = data.R0[0];
     33         IC[1] = data.Rdot0[0];
     34     }
     35 
     36     virtual void compute(const Tinput& U, Trhs& rhs, Real const t, void const* const data=nullptr)
     37     {
     38         const BubbleData& bd = *(BubbleData const* const)data;
     39         const typename Tinput::DataType& u = U[0];
     40         typename Trhs::DataType& r = rhs[0];
     41 
     42         assert(Tinput::DataType::SIZE == 2);
     43         assert(u[0] > 0.0);
     44         assert(data != nullptr);
     45         assert(bd.Nbubbles == 1);
     46 
     47         r[0] = u[1];
     48 
     49         r[1]  = bd.rInv*(bd.pv - bd.p0 - bd.pext(t));
     50         r[1] += bd.rInv*(bd.p0 - bd.pv + static_cast<Real>(2.0)*bd.sigma/bd.R0[0])*std::pow(bd.R0[0]/u[0], 3.0*bd.gamma);
     51         r[1] -= static_cast<Real>(1.5)*u[1]*u[1];
     52         r[1] -= static_cast<Real>(4.0)*bd.nuL*u[1]/u[0];
     53         r[1] -= static_cast<Real>(2.0)*bd.sigma*bd.rInv/u[0];
     54         r[1] /= u[0];
     55 
     56         assert(!isnan(r[0]));
     57         assert(!isnan(r[1]));
     58     }
     59 
     60     virtual void write(const size_t step, const Real t, const Real dt, const Tinput& U, const void * const data=nullptr)
     61     {
     62         const BubbleData& bd = *(BubbleData const* const)data;
     63         const Real odat[4] = {U[0][0], U[0][1], bd.pext(t), bd.pBubble(U)};
     64         m_dumper->write(step, t, dt, &odat[0], 4);
     65         // m_dumper->write(step, t, dt, U.data(), Tinput::DataType::SIZE*U.size());
     66     }
     67 };
     68 
     69 #endif /* RAYLEIGHPLESSET_H_DLCI4XQ3 */