Gabriel WlazłowskiThe st-wslda codes support the following initialization modes:
inittype 0 # 0 - create uniform solution and start from it
# 10 - read uniform solution from file `inprefix`/uniform.solution
#
# 5 - start from st-wslda checkpoint, inprefix points to folder with checkpoint binary file
#
# -1 - custom initialization- Initializing solver with a uniform solution
- Starting from checkpoint
- Custom initialization of the solver
Initializing solver with a uniform solution
The simplest way to start the self-consistent process is to initialize buffers with a solution of a uniform problem. To generate the uniform solution user should use:
inittype 0 # 0 - create uniform solution and start from it Alternatively, if a uniform solution is already existing the user can read it from file:
inittype 10 # 10 - read uniform solution from file `inprefix`/uniform.solution
inprefix test # uniform solution is stored in this location Following flags in the input file are used to control the generation process of the uniform solution:
# ----------- IN CASE OF: inittype=0 ----------------
# init0Na 28.0 # requested number of particles of type a, default=Na
# init0Nb 28.0 # requested number of particles of type b, default=Nb
# init0muchange 0.1 # change rate of chemical potentials, default muchange
# init0Tstart 0.1 # start temperature, in units of eF, default value is equal to init0Tstop
# init0Tstop 0.01 # Temperature, in units of eF, default=temperature
# init0DeltaT 0.01 # change of temperature in units of eF, default 0.01
# init0eps 1.0e-9 # epsilon for convergence, default 1.0e-9
# init0scmix 0.5 # mixing parameter in self-consistent process, default=linearmixing
# init0maxiter 1000 # maximum number of iterations, default=10000
init0debug 1 # debug level, default=0 (no debug info), 1 (basic debug info), 2 (detailed debug info)
init0save 1 # 0 - no saving of solution to file, 1 - save solution to the file Note that, many of them are commented out, and then the values are taken from corresponding tags of the main code.
For spin symmetric systems to code typically converges without any problems. For spin-imbalanced systems, the uniform solution typically is hard to obtain. To overcome this problem, we introduced a method where we start from the solution at finite temperature init0Tstart and we gradually decrease the temperature to destination temperature init0Tstop with step init0DeltaT. Such a problematic case is shown below (for temperature 0.01):
# CREATING UNIFORM SOLUTION...
# DEBUG: n_a=0.017708, n_b=0.018750
# DEBUG: eF_a=0.516086, eF_b=0.536132, eF_avg=0.526157
# DEBUG: N_a=17.000000, N_b=18.000000
# DEBUG: p=-0.028571
# DEBUG: alph_a=1.089111, alph_b=1.098021, alph_plus=1.093566
# DEBUG: dalphm_dna=4.393928, dalphm_dnb=-4.149821
# DEBUG: dalphp_dna=0.856253, dalphp_dnb=-0.808683
# DEBUG: dtildeC_dna=-0.323965, dtildeC_dnb=-0.274273
# DEBUG: dD_dna=-0.289825, dD_dnb=-0.265175
# DEBUG: D=-0.006063, tC=-0.032638
# DEBUG: tau_a=0.010967, tau_b=0.012063
# DEBUG: delta=0.258043, nu=0.007702
# TEMPCONV: T=0.010000, iter=10000, delta/eF_a=0.509012, mu_a/eF_a=-0.061340, delta/eF_b=0.489980, mu_b/eF_b=0.911447
# WARNING: MAXITER REACHED!
# TEMPCONV: T=0.010000, energy_kin=1.923270, energy_pot=-0.526501, energy_pair=-0.967234, energy_tot=0.429535
# ENERGY CUT-OFF: EC1=5.150865
# ENERGY CUT-OFF: EC2=-4.682460
# ENERGY CUT-OFF: EC=4.909640
# UNIFORM SOLUTION: delta/eF_a= 0.5090, mu_a/eF_a= -0.0613, delta/eF_b= 0.4900, mu_b/eF_b= 0.9114, ec= 4.9096
# UNIFORM SOLUTION: energy_kin= 1.923270037681, energy_pot= -0.526501268447, energy_pair= -0.967234017120, energy_tot= 0.429534752115
# UNIFORM SOLUTION: nwf=970It was obtained for default values, and the code indicates # WARNING: MAXITER REACHED! which means that the self-consistency was not obtained. However, when setting:
init0Tstart 0.1 # start temperature, in units of eF, we obtain:
# CREATING UNIFORM SOLUTION...
# DEBUG: n_a=0.017708, n_b=0.018750
# DEBUG: eF_a=0.516086, eF_b=0.536132, eF_avg=0.526157
# DEBUG: N_a=17.000000, N_b=18.000000
# DEBUG: p=-0.028571
# DEBUG: alph_a=1.089111, alph_b=1.098021, alph_plus=1.093566
# DEBUG: dalphm_dna=4.393928, dalphm_dnb=-4.149821
# DEBUG: dalphp_dna=0.856253, dalphp_dnb=-0.808683
# DEBUG: dtildeC_dna=-0.323965, dtildeC_dnb=-0.274273
# DEBUG: dD_dna=-0.289825, dD_dnb=-0.265175
# DEBUG: D=-0.006063, tC=-0.032638
# DEBUG: tau_a=0.010967, tau_b=0.012063
# DEBUG: delta=0.258043, nu=0.007702
# TEMPCONV: T=0.100000, iter=2911, delta/eF_a=0.503474, mu_a/eF_a=0.195268, delta/eF_b=0.484650, mu_b/eF_b=0.668692
# TEMPCONV: T=0.100000, energy_kin=1.915994, energy_pot=-0.526501, energy_pair=-0.948608, energy_tot=0.440885
# TEMPCONV: T=0.090000, iter=1873, delta/eF_a=0.503667, mu_a/eF_a=0.163986, delta/eF_b=0.484835, mu_b/eF_b=0.698980
# TEMPCONV: T=0.090000, energy_kin=1.916098, energy_pot=-0.526501, energy_pair=-0.949330, energy_tot=0.440267
# TEMPCONV: T=0.080000, iter=1482, delta/eF_a=0.503698, mu_a/eF_a=0.132936, delta/eF_b=0.484865, mu_b/eF_b=0.729129
# TEMPCONV: T=0.080000, energy_kin=1.915824, energy_pot=-0.526501, energy_pair=-0.949435, energy_tot=0.439888
# TEMPCONV: T=0.070000, iter=1130, delta/eF_a=0.503659, mu_a/eF_a=0.102453, delta/eF_b=0.484828, mu_b/eF_b=0.758759
# TEMPCONV: T=0.070000, energy_kin=1.915372, energy_pot=-0.526501, energy_pair=-0.949276, energy_tot=0.439595
# TEMPCONV: T=0.060000, iter=822, delta/eF_a=0.503586, mu_a/eF_a=0.072729, delta/eF_b=0.484758, mu_b/eF_b=0.787641
# TEMPCONV: T=0.060000, energy_kin=1.914814, energy_pot=-0.526501, energy_pair=-0.948989, energy_tot=0.439323
# TEMPCONV: T=0.050000, iter=561, delta/eF_a=0.503481, mu_a/eF_a=0.043927, delta/eF_b=0.484656, mu_b/eF_b=0.815579
# TEMPCONV: T=0.050000, energy_kin=1.914134, energy_pot=-0.526501, energy_pair=-0.948580, energy_tot=0.439052
# TEMPCONV: T=0.040000, iter=351, delta/eF_a=0.503324, mu_a/eF_a=0.016246, delta/eF_b=0.484505, mu_b/eF_b=0.842323
# TEMPCONV: T=0.040000, energy_kin=1.913266, energy_pot=-0.526501, energy_pair=-0.947983, energy_tot=0.438782
# TEMPCONV: T=0.030000, iter=223, delta/eF_a=0.503090, mu_a/eF_a=-0.010016, delta/eF_b=0.484280, mu_b/eF_b=0.867514
# TEMPCONV: T=0.030000, energy_kin=1.912119, energy_pot=-0.526501, energy_pair=-0.947102, energy_tot=0.438515
# TEMPCONV: T=0.020000, iter=221, delta/eF_a=0.502800, mu_a/eF_a=-0.034271, delta/eF_b=0.484001, mu_b/eF_b=0.890670
# TEMPCONV: T=0.020000, energy_kin=1.910759, energy_pot=-0.526501, energy_pair=-0.946017, energy_tot=0.438242
# TEMPCONV: T=0.010000, iter=112, delta/eF_a=0.502792, mu_a/eF_a=-0.054949, delta/eF_b=0.483994, mu_b/eF_b=0.911258
# TEMPCONV: T=0.010000, energy_kin=1.910379, energy_pot=-0.526501, energy_pair=-0.945957, energy_tot=0.437921
# ENERGY CUT-OFF: EC1=5.146473
# ENERGY CUT-OFF: EC2=-4.682285
# ENERGY CUT-OFF: EC=4.907524
# UNIFORM SOLUTION: delta/eF_a= 0.5028, mu_a/eF_a= -0.0549, delta/eF_b= 0.4840, mu_b/eF_b= 0.9113, ec= 4.9075
# UNIFORM SOLUTION: energy_kin= 1.910378825876, energy_pot= -0.526501268447, energy_pair= -0.945956681969, energy_tot= 0.437920875460
# UNIFORM SOLUTION: nwf=970Depending on the settings you will need to also modify other init0 parameters.
NOTE: You can check the correctness of settings of the uniform solver by using solve-uniform tool. For more info see here.
Starting from checkpoint
The code can be initialized by an already existing solution generated by st-wslda.
It is the most commonly used option. To activate this option you need to set:
inittype 5 # 5 - start from st-wslda checkpoint, inprefix points to folder with checkpoint binary file
inprefix testUsing checkpoint generated by code of lower dimension
The Toolkit allows for initializing higher dimensional calculations via checkpoint generated by lower-dimensional code, for example, 3D code can be initialized by checkpoints generated by 2D or 1D code.
In the case of resizing the checkpoint file, a message of this type will be provided on stdout:
# CODE: ST-WSLDA-3D
# VERSION: 2021.02.22
# LATTICE: 8 x 10 x 12
...
# INSPECTING CHECKPOINT FILE `test/checkpoint.dat`
# CHECKPOINT FOR 2D LATTICE: [NX,NY,NZ]=[8,10,12], [DX,DY,DZ]=[1.000,1.000,1.000], [LX,LY,LZ]=[8.000,10.000,12.000]
# !!! --- WARNING --- WARNING --- WARNING --- WARNING --- WARNING --- WARNING --- !!!
# Dimensonality of the lattice has changed!
# The code will change dimensionality of given checkpoint data to the new lattice.
# !!! --- ------- --- ------- --- ------- --- ------- --- ------- --- ------- --- !!!
# LOADING CHECKPOINT FILE `test/checkpoint.dat`
# CONVERTING CHECKPOINT FILE: 2D --> 3DUsing checkpoint generated on a lattice with different resolution.
If provided checkpoint file is generated on a lattice with different resolution, than specified in predefines.h, then automatic interpolation will be applied.
Custom initialization of the solver
The static codes can start the self-consistent process from arbitrary initial conditions. To activate this option user must specify in the input file:
inittype -1 # -1 - custom initializationUsing this functionality it is mandatory to initialize value of potentials via modify_potentials(...) routine. Initializtion iteration has index it=-1. Below we present example of such initialization.
void modify_potentials(int it, wslda_density h_densities, wslda_potential h_potentials, double *params, size_t extra_data_size, void *extra_data)
{
// DETERMINE LOCAL SIZES OF ARRAYS (CODE DIMENSIONALITY DEPENDENT)
int lNX=h_densities.nx, lNY=h_densities.ny, lNZ=h_densities.nz; // local sizes
int ix, iy, iz, ixyz;
if(it==-1 && wsldapid==0) printf("SETTING MY VALUES FOR STARTING POINT\n");
// ITERATE OVER ALL POINTS
ixyz=0;
if(it==-1) for(ix=0; ix<lNX; ix++) for(iy=0; iy<lNY; iy++) for(iz=0; iz<lNZ; iz++)
{
double x = DX*(ix-lNX/2);
double y = DY*(iy-lNY/2); // for 1d code y will be always 0
double z = DZ*(iz-lNZ/2); // for 1d and 2d codes z will be always 0
// my custom initialization of the solver
h_densities.rho_a[ixyz]=h_densities.rho_b[ixyz]=0.018; // guess for the density, here I set constant value
// all potentials must be filled, since they are used to construct BdG matrix
h_potentials.delta[ixyz]=0.1 + I*0.0; // guess for delta
h_potentials.alpha_a[ixyz]=h_potentials.alpha_b[ixyz]=1.0; // bare mass
h_potentials.V_a[ixyz]=h_potentials.V_b[ixyz]=0.0; // no mean-field
h_potentials.A_a_x[ixyz]=h_potentials.A_b_x[ixyz]=0.0; // no current potentials
h_potentials.A_a_y[ixyz]=h_potentials.A_b_y[ixyz]=0.0; // no current potentials
h_potentials.A_a_z[ixyz]=h_potentials.A_b_z[ixyz]=0.0; // no current potentials
ixyz++; // go to next point, it should be last line of the triple loop
}
// my guess for chemical potentials
if(it==-1)
{
h_potentials.mu[SPINA]=0.5;
h_potentials.mu[SPINB]=0.5;
}
}