search

SixTrack / sixtracklib

Library for single charged particle simulations in accelerators
GNU Lesser General Public License v2.1
12 stars 16 forks source link

Update delta #139

Open giadarol opened 3 years ago

giadarol commented 3 years ago

Hello, is there a reason why "add_to_energy" affects zeta and "update delta" does not?

Compare: https://github.com/SixTrack/sixtracklib/blob/42b4f06f269602cf464278a81bc6fd589c1d8903/sixtracklib/common/particles.h#L5459

vs

https://github.com/SixTrack/sixtracklib/blob/42b4f06f269602cf464278a81bc6fd589c1d8903/sixtracklib/common/particles.h#L5255

rdemaria commented 3 years ago

Yes, it should update z, maybe it is done through add to energy?

giadarol commented 3 years ago

Yes, it should update z, maybe it is done through add to energy?

It doesn't look like it... Se below:

SIXTRL_INLINE void NS(Particles_update_delta_value)(
    SIXTRL_PARTICLE_ARGPTR_DEC NS(Particles)* SIXTRL_RESTRICT p,
    NS(particle_num_elements_t) const index,
    NS(particle_real_t) const new_delta_value )
{
    typedef NS(particle_real_t) real_t;

    SIXTRL_STATIC_VAR real_t const ONE = ( real_t )1;

    real_t const beta0 = NS(Particles_get_beta0_value)( p, index );
    real_t const delta_beta0 = new_delta_value * beta0;
    real_t const ptau_beta0  = sqrt( delta_beta0 * delta_beta0 +
        ( real_t )2 * delta_beta0 * beta0 + ONE ) - ONE;

    real_t const one_plus_delta = ONE + new_delta_value;
    real_t const rvv    = ( one_plus_delta ) / ( ONE + ptau_beta0 );
    real_t const rpp    = ONE / one_plus_delta;
    real_t const psigma = ptau_beta0 / ( beta0 * beta0 );

    #if !defined( NDEBUG ) && !defined( _GPUCODE )
    SIXTRL_STATIC_VAR real_t const EPS  = ( real_t )1e-9;
    SIXTRL_STATIC_VAR real_t const ZERO = ( real_t )0;

    SIXTRL_ASSERT(   beta0              > ZERO );
    SIXTRL_ASSERT( ( beta0 * beta0    ) > EPS  );
    SIXTRL_ASSERT( ( one_plus_delta   ) > EPS  );
    SIXTRL_ASSERT( ( ONE + ptau_beta0 ) > EPS  );
    SIXTRL_ASSERT( ( delta_beta0 * delta_beta0 +
        ( real_t )2 * delta_beta0 * beta0 + ONE ) > ZERO );

    #endif /* !defined( NDEBUG ) && !defined( _GPUCODE ) */

    NS(Particles_set_delta_value)(  p, index, new_delta_value );
    NS(Particles_set_rvv_value)(    p, index, rvv );
    NS(Particles_set_rpp_value)(    p, index, rpp );
    NS(Particles_set_psigma_value)( p, index, psigma );

    return;
}