NOW COMPILES. Improvements to timestepping.

This commit is contained in:
Miikka Vaisala
2020-09-11 17:03:10 +08:00
parent 66662c1603
commit e9e247c526
2 changed files with 6 additions and 7 deletions

View File

@@ -31,7 +31,7 @@
static AcReal timescale = AcReal(1.0);
AcReal
host_timestep(const AcReal& umax, const AcMeshInfo& mesh_info)
host_timestep(const AcReal& umax, const AcReal& vAmax, const AcMeshInfo& mesh_info)
{
const long double cdt = mesh_info.real_params[AC_cdt];
const long double cdtv = mesh_info.real_params[AC_cdtv];
@@ -40,7 +40,7 @@ host_timestep(const AcReal& umax, const AcMeshInfo& mesh_info)
const long double nu_visc = mesh_info.real_params[AC_nu_visc];
const long double eta = mesh_info.real_params[AC_eta];
const long double chi = 0; // mesh_info.real_params[AC_chi]; // TODO not calculated
const long double gamma = mesh_info.real_params[AC_gamma];
const long double gamma = mesh_info.real_params[AC_gamma]; //TODO this does not make sense here at all.
const long double dsmin = mesh_info.real_params[AC_dsmin];
// Old ones from legacy Astaroth
@@ -49,12 +49,11 @@ host_timestep(const AcReal& umax, const AcMeshInfo& mesh_info)
// New, closer to the actual Courant timestep
// See Pencil Code user manual p. 38 (timestep section)
const long double uu_dt = cdt * dsmin / (fabsl(umax) + sqrtl(cs2_sound + 0.0l));
//const long double uu_dt = cdt * dsmin / (fabsl(umax) + sqrtl(cs2_sound + 0.0l));
const long double uu_dt = cdt * dsmin / (fabsl(umax) + sqrtl(cs2_sound + vAmax*vAmax));
const long double visc_dt = cdtv * dsmin * dsmin /
max(max(nu_visc, eta),
max(gamma, chi)); // + 1; // TODO NOTE: comment the +1 out to
// get scientifically accurate results
// MV: White the +1? It was messing up my computations!
max(gamma, chi));
const long double dt = min(uu_dt, visc_dt);
return AcReal(timescale) * AcReal(dt);

View File

@@ -27,6 +27,6 @@
#pragma once
#include "astaroth.h"
AcReal host_timestep(const AcReal& umax, const AcMeshInfo& mesh_info);
AcReal host_timestep(const AcReal& umax, const AcReal& vAmax, const AcMeshInfo& mesh_info);
void set_timescale(const AcReal scale);