<commit_before>/*
 * SWE_Sphere_TS_ln_settls.cpp
 *
 *  Created on: 24 Sep 2019
 *      Author: <NAME> <<EMAIL>>
 *
 *  Changelog:
 *  	2019-10-24: Partly based on plane version
 */

#include "SWE_Sphere_TS_ln_settls_vd.hpp"



bool SWE_Sphere_TS_ln_settls_vd::implements_timestepping_method(const std::string &i_timestepping_method)
{
	/*
	 * Should contain _exp and _settls
	 */
	return (
		!(i_timestepping_method.find("_exp") != std::string::npos)		&&
		(i_timestepping_method.find("_settls") != std::string::npos)	&&
		(i_timestepping_method.find("_vd") != std::string::npos)		&&
		!(i_timestepping_method.find("_only") != std::string::npos)		&&
		true
	);
}


std::string SWE_Sphere_TS_ln_settls_vd::string_id()
{
	return string_id_storage;
}


void SWE_Sphere_TS_ln_settls_vd::run_timestep(
		SphereData_Spectral &io_phi,	///< prognostic variables
		SphereData_Spectral &io_vrt,	///< prognostic variables
		SphereData_Spectral &io_div,	///< prognostic variables

		double i_fixed_dt,
		double i_simulation_timestamp
)
{
	if (timestepping_order == 2)
	{
		run_timestep_2nd_order(io_phi, io_vrt, io_div, i_fixed_dt, i_simulation_timestamp);
		return;
	}
}




void SWE_Sphere_TS_ln_settls_vd::run_timestep_2nd_order(
	SphereData_Spectral &io_U_phi,		///< prognostic variables
	SphereData_Spectral &io_U_vrt,		///< prognostic variables
	SphereData_Spectral &io_U_div,		///< prognostic variables

	double i_dt,		
	double i_simulation_timestamp
)
{
	const SphereData_Spectral &U_phi = io_U_phi;
	const SphereData_Spectral &U_vrt = io_U_vrt;
	const SphereData_Spectral &U_div = io_U_div;

	if (i_simulation_timestamp == 0)
	{
		/*
		 * First time step:
		 * Simply backup existing fields for multi-step parts of this algorithm.
		 */
		U_phi_prev = U_phi;
		U_vrt_prev = U_vrt;
		U_div_prev = U_div;
	}


	/*
	 * Step 1) SL
	 * Compute Lagrangian trajectories based on SETTLS.
	 * This depends on V(t-\Delta t) and V(t).
	 *
	 * See Hortal's paper for equation.
	 */
	SphereData_Physical U_u_lon_prev, U_v_lat_prev;
	ops.vrtdiv_to_uv(U_vrt_prev, U_div_prev, U_u_lon_prev, U_v_lat_prev);

	SphereData_Physical U_u_lon, U_v_lat;
	ops.vrtdiv_to_uv(U_vrt, U_div, U_u_lon, U_v_lat);

	double dt_div_radius = i_dt / simVars.sim.sphere_radius;

	// Calculate departure points
	ScalarDataArray pos_lon_d, pos_lat_d;
	semiLagrangian.semi_lag_departure_points_settls_specialized(
			dt_div_radius*U_u_lon_prev, dt_div_radius*U_v_lat_prev,
			dt_div_radius*U_u_lon, dt_div_radius*U_v_lat,
			pos_lon_d, pos_lat_d		// OUTPUT
	);


	/*
	 * Step 2) Midpoint rule
	 * Put everything together with midpoint rule and solve resulting Helmholtz problem
	 */

	/*
	 * Step 2a) Compute RHS
	 * R = X_D + 1/2 dt L_D + dt N*
	 */

	/*
	 * Compute X_D
	 */
	SphereData_Spectral U_phi_D, U_vrt_D, U_div_D;
	semiLagrangian.apply_sl_timeintegration_vd(
			ops,
			U_phi, U_vrt, U_div,
			pos_lon_d, pos_lat_d,
			U_phi_D, U_vrt_D, U_div_D
		);

	/*
	 * Compute L_D
	 */
	SphereData_Spectral L_U_phi_D, L_U_vrt_D, L_U_div_D;
	if (original_linear_operator_sl_treatment)
	{
		/*
		 * Method 1) First evaluate L, then sample result at departure point
		 */
		const SphereData_Config *sphereDataConfig = io_U_phi.sphereDataConfig;
		SphereData_Spectral L_U_phi(sphereDataConfig, 0), L_U_vrt(sphereDataConfig, 0), L_U_div(sphereDataConfig, 0);

		/*
		 * L_g(U): Linear gravity modes
		 */
		swe_sphere_ts_ln_erk_split_vd->euler_timestep_update_lg(
				U_phi, U_vrt, U_div,
				L_U_phi, L_U_vrt, L_U_div,
				i_simulation_timestamp
			);

		if (coriolis_treatment == CORIOLIS_LINEAR)
		{
			/*
			 * L_c(U): Linear Coriolis effect
			 */
			swe_sphere_ts_ln_erk_split_vd->euler_timestep_update_lc(
					U_phi, U_vrt, U_div,
					L_U_phi, L_U_vrt, L_U_div,
					i_simulation_timestamp
				);
		}

		if (coriolis_treatment == CORIOLIS_SEMILAGRANGIAN)
		{
			semiLagrangian.apply_sl_timeintegration_vd(
					ops,
					L_U_phi, L_U_vrt + coriolis_spectral, L_U_div,		// Coriolis added to vorticity!
					pos_lon_d, pos_lat_d,
					L_U_phi_D, L_U_vrt_D, L_U_div_D
				);

			L_U_vrt_D -= coriolis_spectral;
		}
		else
		{
			semiLagrangian.apply_sl_timeintegration_vd(
					ops,
					L_U_phi, L_U_vrt, L_U_div,
					pos_lon_d, pos_lat_d,
					L_U_phi_D, L_U_vrt_D, L_U_div_D
				);
		}
	}
	else
	{
		/*
		 * Method 2) First get variables on departure points, then evaluate L
		 */

		SphereData_Spectral U_phi_D, U_vrt_D, U_div_D;

		if (coriolis_treatment == CORIOLIS_SEMILAGRANGIAN)
		{
			semiLagrangian.apply_sl_timeintegration_vd(
					ops,
					U_phi, U_vrt + coriolis_spectral, U_div,	// Coriolis added to vorticity!
					pos_lon_d, pos_lat_d,
					U_phi_D, U_vrt_D, U_div_D
				);

			U_vrt_D -= coriolis_spectral;
		}<commit_msg>swapped sign for LC SL<commit_after>/*
 * SWE_Sphere_TS_ln_settls.cpp
 *
 *  Created on: 24 Sep 2019
 *      Author: <NAME> <<EMAIL>>
 *
 *  Changelog:
 *  	2019-10-24: Partly based on plane version
 */

#include "SWE_Sphere_TS_ln_settls_vd.hpp"



bool SWE_Sphere_TS_ln_settls_vd::implements_timestepping_method(const std::string &i_timestepping_method)
{
	/*
	 * Should contain _exp and _settls
	 */
	return (
		!(i_timestepping_method.find("_exp") != std::string::npos)		&&
		(i_timestepping_method.find("_settls") != std::string::npos)	&&
		(i_timestepping_method.find("_vd") != std::string::npos)		&&
		!(i_timestepping_method.find("_only") != std::string::npos)		&&
		true
	);
}


std::string SWE_Sphere_TS_ln_settls_vd::string_id()
{
	return string_id_storage;
}


void SWE_Sphere_TS_ln_settls_vd::run_timestep(
		SphereData_Spectral &io_phi,	///< prognostic variables
		SphereData_Spectral &io_vrt,	///< prognostic variables
		SphereData_Spectral &io_div,	///< prognostic variables

		double i_fixed_dt,
		double i_simulation_timestamp
)
{
	if (timestepping_order == 2)
	{
		run_timestep_2nd_order(io_phi, io_vrt, io_div, i_fixed_dt, i_simulation_timestamp);
		return;
	}
}




void SWE_Sphere_TS_ln_settls_vd::run_timestep_2nd_order(
	SphereData_Spectral &io_U_phi,		///< prognostic variables
	SphereData_Spectral &io_U_vrt,		///< prognostic variables
	SphereData_Spectral &io_U_div,		///< prognostic variables

	double i_dt,		
	double i_simulation_timestamp
)
{
	const SphereData_Spectral &U_phi = io_U_phi;
	const SphereData_Spectral &U_vrt = io_U_vrt;
	const SphereData_Spectral &U_div = io_U_div;

	if (i_simulation_timestamp == 0)
	{
		/*
		 * First time step:
		 * Simply backup existing fields for multi-step parts of this algorithm.
		 */
		U_phi_prev = U_phi;
		U_vrt_prev = U_vrt;
		U_div_prev = U_div;
	}


	/*
	 * Step 1) SL
	 * Compute Lagrangian trajectories based on SETTLS.
	 * This depends on V(t-\Delta t) and V(t).
	 *
	 * See Hortal's paper for equation.
	 */
	SphereData_Physical U_u_lon_prev, U_v_lat_prev;
	ops.vrtdiv_to_uv(U_vrt_prev, U_div_prev, U_u_lon_prev, U_v_lat_prev);

	SphereData_Physical U_u_lon, U_v_lat;
	ops.vrtdiv_to_uv(U_vrt, U_div, U_u_lon, U_v_lat);

	double dt_div_radius = i_dt / simVars.sim.sphere_radius;

	// Calculate departure points
	ScalarDataArray pos_lon_d, pos_lat_d;
	semiLagrangian.semi_lag_departure_points_settls_specialized(
			dt_div_radius*U_u_lon_prev, dt_div_radius*U_v_lat_prev,
			dt_div_radius*U_u_lon, dt_div_radius*U_v_lat,
			pos_lon_d, pos_lat_d		// OUTPUT
	);


	/*
	 * Step 2) Midpoint rule
	 * Put everything together with midpoint rule and solve resulting Helmholtz problem
	 */

	/*
	 * Step 2a) Compute RHS
	 * R = X_D + 1/2 dt L_D + dt N*
	 */

	/*
	 * Compute X_D
	 */
	SphereData_Spectral U_phi_D, U_vrt_D, U_div_D;
	semiLagrangian.apply_sl_timeintegration_vd(
			ops,
			U_phi, U_vrt, U_div,
			pos_lon_d, pos_lat_d,
			U_phi_D, U_vrt_D, U_div_D
		);

	/*
	 * Compute L_D
	 */
	SphereData_Spectral L_U_phi_D, L_U_vrt_D, L_U_div_D;
	if (original_linear_operator_sl_treatment)
	{
		/*
		 * Method 1) First evaluate L, then sample result at departure point
		 */
		const SphereData_Config *sphereDataConfig = io_U_phi.sphereDataConfig;
		SphereData_Spectral L_U_phi(sphereDataConfig, 0), L_U_vrt(sphereDataConfig, 0), L_U_div(sphereDataConfig, 0);

		/*
		 * L_g(U): Linear gravity modes
		 */
		swe_sphere_ts_ln_erk_split_vd->euler_timestep_update_lg(
				U_phi, U_vrt, U_div,
				L_U_phi, L_U_vrt, L_U_div,
				i_simulation_timestamp
			);

		if (coriolis_treatment == CORIOLIS_LINEAR)
		{
			/*
			 * L_c(U): Linear Coriolis effect
			 */
			swe_sphere_ts_ln_erk_split_vd->euler_timestep_update_lc(
					U_phi, U_vrt, U_div,
					L_U_phi, L_U_vrt, L_U_div,
					i_simulation_timestamp
				);
		}

		if (coriolis_treatment == CORIOLIS_SEMILAGRANGIAN)
		{
			semiLagrangian.apply_sl_timeintegration_vd(
					ops,
					L_U_phi, L_U_vrt - coriolis_spectral, L_U_div,		// Coriolis added to vorticity!
					pos_lon_d, pos_lat_d,
					L_U_phi_D, L_U_vrt_D, L_U_div_D
				);

			L_U_vrt_D += coriolis_spectral;
		}
		else
		{
			semiLagrangian.apply_sl_timeintegration_vd(
					ops,
					L_U_phi, L_U_vrt, L_U_div,
					pos_lon_d, pos_lat_d,
					L_U_phi_D, L_U_vrt_D, L_U_div_D
				);
		}
	}
	else
	{
		/*
		 * Method 2) First get variables on departure points, then evaluate L
		 */

		SphereData_Spectral U_phi_D, U_vrt_D, U_div_D;

		if (coriolis_treatment == CORIOLIS_SEMILAGRANGIAN)
		{
			semiLagrangian.apply_sl_timeintegration_vd(
					ops,
					U_phi, U_vrt - coriolis_spectral, U_div,	// Coriolis added to vorticity!
					pos_lon_d, pos_lat_d,
					U_phi_D, U_vrt_D, U_div_D
				);

			U_vrt_D += coriolis_spectral;
		}
{'fraction_non_alphanumeric': 0.09924556370205079, 'fraction_numerical': 0.005950483476782489, 'mean_word_length': 4.046648793565684, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6542423', 'n_tokens_mistral': 4683, 'n_tokens_neox': 4012, 'n_words': 684}