Skip to content

GitLab

  • Menu
Projects Groups Snippets
    • Loading...
  • Help
    • Help
    • Support
    • Community forum
    • Submit feedback
    • Contribute to GitLab
  • Sign in
  • wslda wslda
  • Project information
    • Project information
    • Activity
    • Labels
    • Members
  • Repository
    • Repository
    • Files
    • Commits
    • Branches
    • Tags
    • Contributors
    • Graph
    • Compare
  • Issues 0
    • Issues 0
    • List
    • Boards
    • Service Desk
    • Milestones
  • Merge requests 0
    • Merge requests 0
  • CI/CD
    • CI/CD
    • Pipelines
    • Jobs
    • Schedules
  • Deployments
    • Deployments
    • Environments
    • Releases
  • Monitor
    • Monitor
    • Incidents
  • Packages & Registries
    • Packages & Registries
    • Package Registry
    • Container Registry
    • Infrastructure Registry
  • Analytics
    • Analytics
    • CI/CD
    • Repository
    • Value stream
  • Wiki
    • Wiki
  • Snippets
    • Snippets
  • Activity
  • Graph
  • Create a new issue
  • Jobs
  • Commits
  • Issue Boards
Collapse sidebar
  • wtools
  • wsldawslda
  • Wiki
  • Tracking of selected states

Last edited by Gabriel Wlazłowski Feb 26, 2024
Page history

Tracking of selected states

  • Tracking the states using build in mechanism
    • Defining the subset
    • Limitations
    • Example
  • Writing to files selected states
  • Tracking the state using write_wave_functions(...)

Tracking the states using build in mechanism

VERSION>=2021.09.18

Defining the subset

W-SLDA Toolkit allows the tracking of observables arising from selected states in time. Presently, the following quantities can be tracked:

  • subset_rho_a:
    n_{\uparrow}^{\textrm{(subset)}}(r)= \sum_{E_{\textrm{min}}^{\textrm{(subset)}}<E_n<E_{\textrm{max}}^{\textrm{(subset)}}}|u_{n,\uparrow}(r)|^2 f_{\beta}(E_n)
  • subset_rho_b:
    n_{\downarrow}^{\textrm{(subset)}}(r) = \sum_{E_{\textrm{min}}^{\textrm{(subset)}}<E_n<E_{\textrm{max}}^{\textrm{(subset)}}}|v_{n,\downarrow}(r)|^2 f_{\beta}(-E_n)
  • subset_j_a_x, subset_j_a_y, subset_j_a_z:
    \vec{j}_{\uparrow}^{\textrm{(subset)}}(r) = -\sum_{E_{\textrm{min}}^{\textrm{(subset)}}<E_n<E_{\textrm{max}}^{\textrm{(subset)}}} \textrm{Im}[u_{n,\uparrow}(r)\nabla u_{n,\uparrow}^*(r)]f_{\beta}(E_n)
  • subset_j_b_x, subset_j_b_y, subset_j_b_z:
    \vec{j}_{\downarrow}^{\textrm{(subset)}}(r) = \sum_{E_{\textrm{min}}^{\textrm{(subset)}}<E_n<E_{\textrm{max}}^{\textrm{(subset)}}} \textrm{Im}[v_{n,\downarrow}(r)\nabla v_{n,\downarrow}^*(r)]f_{\beta}(-E_n)

To control interval range E_{\textrm{min}}^{\textrm{(subset)}}<E_n<E_{\textrm{max}}^{\textrm{(subset)}} use input file variables:

subsetMinEn             0.0     # in eF units, deafault=0.0
subsetMaxEn             0.0     # in eF units, deafault=0.0

In the case of spin-imbalanced systems, it is also convenient to introduce a shift of quasiparticle energies when selecting the states: E_n \rightarrow E_n+\frac{\mu_{\uparrow}-\mu_{\downarrow}}{2}. This can be done automatically by enabling

subsetShiftDmu          1       # if 1 then apply extra shift of quasiparticle energies by (mu_a-mu_b)/2, default=0

Limitations

The states to be tracked and contracted into subset_rho and subset_j are tagged at the beginning of the simulations. We select states based on their energies at t=0. However, as the dynamics proceed, the energies of the states can change as well, so some of them can acquire energy that is out of the given range, and some that were not included initially can acquire energy within the considered energy range. These effects are not taken into account by the presented prescription.

Example

Below we present a snapshot from a simulation with SLDA functional for the spin-symmetric system. td-wslda-2d has been used. The unitary Fermi gas is confined in a tube, and the vortex dipole is imprinted in the initial state. Following options for subset tracking were used

# ---------------- SUBSET TRACKING ------------------
# See: Wiki -> Tracking of selected states
# observables (densities, currents) arising from state in the energy interval En in [subsetMinEn,subsetMaxEn]
# will be computed, if  subsetMinEn=subsetMaxEn=0  the functionality is disables
subsetMinEn             -0.35     # in eF units, deafault=0.0
subsetMaxEn             0.35     # in eF units, deafault=0.0
subsetShiftDmu          1       # if 1 then apply extra shift of quasiparticle energies by (mu_a-mu_b)/2, default=0

For the unitary Fermi gas \Delta/\varepsilon_{F}=0.5, thus, as the subset, we select only in-gap states, which in this case corresponds to Andreev vortex states. They are localized in the vortex core. It is visualized below: the left half shows \Delta(r)/\varepsilon_{F}, while the right half shows subset_rho_a. subset

Writing to files selected states

VERSION>=2024.01.31

Tracking the state using write_wave_functions(...)

The use can track the selected state via write_wave_functions(...) function located in logger.h. The example below shows an example code that writes periodically to files states with quasiparticle energies smaller than the given threshold value. The example is for 2D code.

static int wffileid; // line id
static int wfcall;   // call id of the function
int *idx_wf_for_storage;
static int cnt_wf_for_storage;
int write_wave_functions(int it, int nxyz, int nwfip, int nwf, double beta, MPI_Comm comm,
                         double complex *h_wf, double *h_qpe, double *h_kky, double *h_kkz, int *h_cnt,  // Host pointers
                         double complex *d_wf,                                                           // Device pointers
                         double *params, size_t extra_data_size, void *extra_data
                        )
{
    // settings mapping
    int sweep = params[11];
    double delta = params[12];

    // do periodically
    wfcall++;
    if((wfcall-1) % sweep != 0) return 0; // break

    // example code demonstrating writing wave functions to data format (readable by VisIt).
    int ip, np;
    MPI_Comm_size(comm, &np); // np = total number of processes
    MPI_Comm_rank(comm, &ip); // id of process st 0 <= ip < np

    // printf("S--> %d %d %d %d\n", ip, np, nwfip, nwf);

    // Copy wave functions from device to host
    memcopy_gpu2host(d_wf, h_wf, (size_t)2*nxyz*nwfip*sizeof(double complex));
    // set pointers to u and v components
    double complex *wf_u=h_wf+(size_t)0*nxyz*nwfip; // pointer to nwfip wave functions
    double complex *wf_v=h_wf+(size_t)1*nxyz*nwfip; // pointer to nwfip wave functions

    int i, ixy, cnt ;

    // identify only states that have e_n<Delta, where Delta is stored in delta
    // and save their indices to the global array
    // do it only for the firts call of the function
    if(wffileid==0)
    {
        cnt=0;
        for(i=0; i<nwfip; i++) if(fabs(h_qpe[i])<delta) {cnt++;}
        cnt_wf_for_storage = cnt;

        if(cnt_wf_for_storage>0) idx_wf_for_storage = (int *)malloc(cnt_wf_for_storage*sizeof(int));
        else idx_wf_for_storage=NULL;

        cnt=0;
        for(i=0; i<nwfip; i++) if(fabs(h_qpe[i])<delta) {idx_wf_for_storage[cnt]=i; cnt++;}
    }

    double *my_qpe=NULL, *my_kkz=NULL;
    if(cnt_wf_for_storage>0)
    {
        my_qpe = (double *)malloc(sizeof(double)*cnt_wf_for_storage);
        my_kkz = (double *)malloc(sizeof(double)*cnt_wf_for_storage);
    }

    for(i=0; i<cnt_wf_for_storage; i++) my_qpe[i]=h_qpe[idx_wf_for_storage[i]]; // make local copy
    for(i=0; i<cnt_wf_for_storage; i++) my_kkz[i]=h_kkz[idx_wf_for_storage[i]]; // make local copy
    for(i=0; i<cnt_wf_for_storage; i++) for(ixy=0; ixy<NXY; ixy++) wf_u[i*nxyz+ixy] = wf_u[idx_wf_for_storage[i]*nxyz+ixy]; // reorganize
    for(i=0; i<cnt_wf_for_storage; i++) for(ixy=0; ixy<NXY; ixy++) wf_v[i*nxyz+ixy] = wf_v[idx_wf_for_storage[i]*nxyz+ixy]; // reorganize
    nwfip=cnt_wf_for_storage;
    MPI_Allreduce(&nwfip, &nwf, 1, MPI_INT, MPI_SUM, comm);
    if(ip==0) printf("# [write_wave_functions] IDENTIFIED %d STATES WITH |E_n|<%f\n", nwf, delta);

    // get info about the number of wave functions kept by each MPI process
    int * wf_tbl = (int *)malloc(np*sizeof(int));
    MPI_Gather( &nwfip , 1 , MPI_INT , wf_tbl , 1 , MPI_INT , 0 , MPI_COMM_WORLD ) ;
    MPI_Bcast( wf_tbl , np , MPI_INT , 0 , MPI_COMM_WORLD ) ;

    // use MPI I/O to write a file
    char file_name[128];
    MPI_File fh;
    MPI_Offset my_offset=0;
    MPI_Status status;

    for(i=0; i<ip; i++) my_offset+=sizeof(double complex)*nxyz*wf_tbl[i];

    // u component
    sprintf(file_name, "%s_%06d_wf_u.wdat", md.outprefix, wffileid);
    if(ip==0) printf("# [write_wave_functions] WRITING U COMPONENTS TO FILE `%s`\n", file_name);
    MPI_File_open(comm, file_name, MPI_MODE_CREATE | MPI_MODE_WRONLY, MPI_INFO_NULL, &fh);
    MPI_File_seek(fh, my_offset, MPI_SEEK_SET);
    MPI_File_write_all(fh, wf_u , nxyz*nwfip, MPI_DOUBLE_COMPLEX, &status);
    MPI_Get_count(&status, MPI_DOUBLE_COMPLEX, &cnt);
    if(cnt!=nxyz*nwfip) { printf("# PROBLEM WITH WRITING OF U COMPONENTS BY PROCESS %s!\n", ip); return 1; }
    MPI_File_close(&fh);

    // v component
    sprintf(file_name, "%s_%06d_wf_v.wdat", md.outprefix, wffileid);
    if(ip==0) printf("# [write_wave_functions] WRITING V COMPONENTS TO FILE `%s`\n", file_name);
    MPI_File_open(comm, file_name, MPI_MODE_CREATE | MPI_MODE_WRONLY, MPI_INFO_NULL, &fh);
    MPI_File_seek(fh, my_offset, MPI_SEEK_SET);
    MPI_File_write_all(fh, wf_v , nxyz*nwfip, MPI_DOUBLE_COMPLEX, &status);
    MPI_Get_count(&status, MPI_DOUBLE_COMPLEX, &cnt);
    if(cnt!=nxyz*nwfip) { printf("# PROBLEM WITH WRITING OF V COMPONENTS BY PROCESS %s!\n", ip); return 2; }
    MPI_File_close(&fh);

    // to be able to see the wave functions in VisIt add wtxt file
    // see documentation of WDATA for more details
    wdata_metadata wdmd;
    wdmd.datadim = CODEDIM;
    wdmd.nx = NX; wdmd.ny = NY; wdmd.nz = NZ;
    wdmd.dx = DX; wdmd.dy = DY; wdmd.dz = DZ;
    sprintf(wdmd.prefix, "%s_%06d", md.outprefix, wffileid);
    wdmd.t0 = 0.0; wdmd.dt = 1.0; wdmd.cycles=nwf; wdmd.nvar=0;
    wdata_variable vwf_u = {"wf_u", "complex", "none", "wdat"}; wdata_add_variable(&wdmd, &vwf_u);
    wdata_variable vwf_v = {"wf_v", "complex", "none", "wdat"}; wdata_add_variable(&wdmd, &vwf_v);
    sprintf(file_name, "%s_%06d.wtxt", md.outprefix, wffileid);
    if(ip==0) wdata_write_metadata_to_file(&wdmd,file_name);

    // write En and kz
    my_offset=0;
#define __LINE_LGHT 128
    char line[__LINE_LGHT];
    sprintf(line, "%5d %12.8f %12.8f \n", 0, 0.0, 0.0); // line format
    int _line_lgh = strlen(line); // test its lenght to get info about number of bytes per line
    for(i=0; i<ip; i++) my_offset+=sizeof(char)*_line_lgh*wf_tbl[i];
    int startidx=0;
    for(i=0; i<ip; i++) startidx+=wf_tbl[i];

    sprintf(file_name, "%s_%06d_Ekz.txt", md.outprefix, wffileid);
    if(ip==0) printf("# [write_wave_functions] WRITING En AND kz TO FILE `%s`\n", file_name);
    MPI_File_open(comm, file_name, MPI_MODE_CREATE | MPI_MODE_WRONLY, MPI_INFO_NULL, &fh);
    for(i=0; i<nwfip; i++)
    {
        sprintf(line, "%5d %12.8f %12.8f \n", startidx+i, my_kkz[i], my_qpe[i]);  // prepare line
        MPI_File_write_at(fh, my_offset, line , _line_lgh, MPI_CHAR, &status);     // write it
        MPI_Get_count(&status, MPI_CHAR, &cnt);                                    // get status
        if(cnt!=_line_lgh) { printf("# PROBLEM WITH WRITING OF En AND kz BY PROCESS %s!\n", ip); return 3; }
        my_offset+=sizeof(char)*_line_lgh;                                         // move pointer to next line
    }
    MPI_File_close(&fh);
#undef __LINE_LGHT

    // clear memory
    free(wf_tbl);

    wffileid++; // store next id
    free(my_kkz);
    free(my_qpe);

    return 0; // return OK status
}
Clone repository

Content of Documentation
Official webpage
W-BSK Toolkit