//  DISP_LAT.C
//  © 2021 Peter J. Meyer

#include <process.h>
#include <math.h>

#include "iss.h"

#define SITE_IN_SPANNING_CLUSTER_2(i,j)     ( sites2[i][j] & IN_SPANNING_CLUSTER )
#define SITE_IN_SPANNING_CLUSTER_3(i,j,k)   ( sites2[i][j][k] & IN_SPANNING_CLUSTER )
#define SITE_IN_SPANNING_CLUSTER_4(i,j,k,l) ( sites2[i][j][k][l] & IN_SPANNING_CLUSTER )

#define ENABLE_CONTRACTED_MAP true
#define MAX_MAP_DIMENSIONALITY 2

#define I ((3*i)+1)
#undef  J
#define J ((3*j)+1)

extern double v_bond_probability;

extern short int w_i, w_j, w_k, w_l;

int print_contracted_map = ENABLE_CONTRACTED_MAP;

static int max_size_for_map[] = { 0, 0, 24, 12, 6 };
//  Maximum sizes for maps based on dimensionality.

static char slice[72][78];

static int q_state_potts_spins[MAX_Q_VALUE+1];
static int total_q_state_potts_spins;

static void blank_map_slice(void);
static void zero_potts_spin_counts(void);
static void write_map_slice(FILE *file, int context, int v_bonds);
static void write_2d_map(FILE *file, int context, int v_bonds);
static void write_3d_map(FILE *file, int context, int v_bonds);
static void write_4d_map(FILE *file, int context, int v_bonds);

/*----------------------------*/
void display_map(char *filepath,
                 int context)       //  0 = general purpose
                                    //  1 = display_spanning_path
                                    //  2 = Swendsen-Wang before cluster flip
                                    //  3 = Swendsen-Wang after cluster flip
                                    //  4 = Wolff after cluster flip
{
FILE *map_file;

if ( !open_file(filepath,"wt",&map_file) )
    display_file_error_message("open","map",filepath);  //  Does not return.

if ( dimensionality > MAX_MAP_DIMENSIONALITY )
    fprintf(map_file,"\nMap file not currently supported for d>%d.\n",
    MAX_MAP_DIMENSIONALITY);
else if ( size > max_size_for_map[dimensionality] )
    fprintf(map_file,"\nMaximum size for map of %dd lattice is %d.\n",
        dimensionality,max_size_for_map[dimensionality]);
else
    {
    switch ( dimensionality )
        {
        case 2:
        //  Write bonds map.
        write_2d_map(map_file,0,false);
        if ( context >= 2 )    
            //  Write virtual bonds map.
            write_2d_map(map_file,context,true);
        break;

        case 3:
        //  write_3d_map(map_file,write_model_info);
        break;

        case 4:
        //  write_4d_map(map_file,write_model_info);
        break;
        }
    }

fclose(map_file);

//  Call up NotePad (rather than WordPad) to display the map file.
spawnlp(P_NOWAIT,"notepad.exe","notepad.exe",filepath,NULL);

//  printf("\nPress a key ... ");
if ( getch() == ESCAPE )
    exit(0);
}


/*--------------------------------*/
static void write_2d_map(FILE *file,
                         int context,
                         int v_bonds)
{
int i, j, occupied, spin;
unsigned char bonds;

blank_map_slice();
zero_potts_spin_counts();

for ( i=0; i<size; i++ )
    {
    for ( j=0; j<size; j++ )
        {
        occupied = spin = spin_at_site(i,j);
        if ( occupied )
            {
            if ( v_bonds )
                bonds = v_bonds2[i][j];
            else
                bonds = bonds2[i][j];
            if ( !honeycomb_lattice )
                {
                if ( bonds & BIT1_0 )
                    slice[I+1][J] = '|';
                    //  slice[I+1][J] = '0' + (cluster_numbers2[i][j]%10);
                if ( bonds & BIT1_1 )
                    slice[I-1][J] = '|';
                if ( bonds & BIT1_2 )
                    slice[I][J+1] = '-';
                if ( bonds & BIT1_3 )
                    slice[I][J-1] = '-';
                if ( bonds & BIT1_4 )
                    slice[I+1][J+1] = '\\';
                if ( bonds & BIT1_5 )
                    {
                    slice[I-1][J-1] = '\\';
                    //  Don't print contracted map since
                    //  this would not show these characters.
                    print_contracted_map = false;
                    }
                if ( bonds & BIT1_6 )
                    slice[I+1][J-1] = '/';
                if ( bonds & BIT1_7 )
                    {
                    slice[I-1][J+1] = '/';
                    print_contracted_map = false;
                    }
                }
            else
                {
                if ( major_axis == 0 )
                    {
                    if ( bonds & BIT1_0 )
                        slice[I+1][J] = '|';
                    if ( bonds & BIT1_1 )
                        slice[I-1][J] = '|';
                    if ( bonds & BIT1_2 )
                        {
                        if ( (i+j)%2 )
                            slice[I][J-1] = '-';
                        else
                            slice[I][J+1] = '-';
                        }
                    }
                else  //  major_axis = 1
                    {
                    if ( bonds & BIT1_0 )
                        slice[I][J+1] = '-';
                    if ( bonds & BIT1_1 )
                        slice[I][J-1] = '-';
                    if ( bonds & BIT1_2 )
                        {
                        if ( (i+j)%2 )
                            slice[I-1][J] = '|';
                        else
                            slice[I+1][J] = '|';
                        }
                    }
                }
            if ( goal_is_percolation_threshold )
                {
                if ( SITE_IN_SPANNING_CLUSTER_2(i,j) )
                    slice[I][J] = 'S';
                else if ( SITE_IS_VISITED_2(i,j) )
                    slice[I][J] = 'V';
                else
                    slice[I][J] = 'O';
                }
            else if ( model_is_ising )
                slice[I][J] = ( spin == UP_SPIN ? 'X' : '0' );
            else  //  model is q-state Potts
                {
                q_state_potts_spins[spin]++;
                total_q_state_potts_spins++;
                if ( q < 10 )
                    slice[I][J] = spin + '0';
                else
                    {
                    slice[I][J] = spin%10 + '0';
                    if ( spin > 9 )
                        slice[I][J-1] = spin/10 + '0';
                    }
                }
            }
        }
    }

write_map_slice(file, context, v_bonds);
}

/*------------------------------------*/
static void zero_potts_spin_counts(void)
{
int i;

for ( i=1; i<=q; i++ )
    q_state_potts_spins[i] =0;

total_q_state_potts_spins = 0;    
}

/*-----------------------------*/
static void blank_map_slice(void)
{
int i, j;

for ( i=0; i<3*size; i++ )
    for ( j=0; j<3*size; j++ )
        slice[i][j] = SPACE;

//  Put dots around the boundaries.
for ( j=0; j<3*size; j++ )
    slice[0][j] = '.';
for ( j=1; j<3*size-1; j++ )
    slice[3*size-1][j] = '.';
for ( i=1; i<3*size; i++ )
    slice[i][0] = slice[i][3*size-1] = ':';
}

/*-----------------------------------*/
static void write_map_slice(FILE *file, 
                            int context,
                            int v_bonds)
{
int i, j;
double magn, site_pc, bond_pc;

if ( context < 2 )
    {

    fprintf(file,"\nInput file: %s.  Runtime: %s %s",
        input_filepath,system_date(temp,0),system_time(temp2,0));

    fprintf(file,"\n\nlattice dimensionality/type/size:  %dd  %s (%d)  %d",
        dimensionality,lattice_types[l_type].name,
        lattice_types[l_type].coord_num,size);
    for ( i=2; i<=dimensionality; i++ )
        fprintf(file,"x%d",size);
    if ( direction_table_is_parity_dependent && MAJOR_AXIS_SPECIFICATION_PERMITTED )
        fprintf(file,"\nmajor axis:          %d",major_axis);

    #if FREE_BOUNDARY_CONDITIONS_POSSIBLE
    fprintf(file,"\n%s boundary conditions",
        ( free_boundary_conditions ? "free" : "periodic" ));
    #endif

    if ( seed )
        fprintf(file,"    Random seed = %u",seed);

    fprintf(file,"\ndynamics:           %s",dynamics_types[d_type][1]);
    if ( dynamics_is_swendsen_wang || dynamics_is_wolff )
        fprintf(file," (p = %.5f)",v_bond_probability);
    else
        fprintf(file,"\nspin selection:     %s",spin_selns[ss_type][1]);

    fprintf(file,"\n\nsite concentration = %.6f, bond concentration = %.6f",
        site_concentration, bond_concentration);

    site_pc = ((double)num_spins*100)/num_sites;

    fprintf(file,"\n%s of %s lattice sites occupied (%.2f%%).\n",
        ultoa_commas(num_spins,temp),
        ultoa_commas(num_sites,temp1),site_pc);

    bond_pc = ((double)num_opened_bonds*100)/num_pairs_nn;

    fprintf(file,"%s bonds opened (%.2f%% of the %s pairs of nn sites).\n",
            ultoa_commas(num_opened_bonds,temp),bond_pc,ultoa_commas(num_pairs_nn,temp2));
    }

if ( goal_is_percolation_threshold )
    fprintf(file,"\nS = site on spanning path, "
                   "V = other visited site, O = other occupied site");
else
    {
    if ( context < 2 )
        {
        if ( model_is_ising )
            fprintf(file,"\nIsing model, ");
        else 
            fprintf(file,"\n%d-state Potts model, ",q);
        fprintf(file,"%dd %s lattice",
            dimensionality,lattice_types[l_type].name);
        if ( direction_table_is_parity_dependent && MAJOR_AXIS_SPECIFICATION_PERMITTED )
            fprintf(file,", major axis = %d",major_axis);
        fprintf(file,"\ninitial magnetization = %.3f, temperature = %.3f",
            initial_magnetization, temperature);
        magn = magnetization_per_spin();
        fprintf(file,"\nMcs = %-6lu ",mcs_num);

        if ( !num_spin_flips_overflow )
            fprintf(file,"num_spin_flips = %s",
                ultoa_commas((unsigned long)num_spin_flips,temp));
        else
            fprintf(file,"Number of millions of spin flips = %.0f",
                (num_spin_flips + pow(2,32)*num_spin_flips_overflow)/1e6);

        fprintf(file,"\nAbs.mag.= %-6d   Mag.=%6.3f  ",abs_magnetization,magn);
        if ( internal_energy_measured )
            fprintf(file,"Abs.int.en.= %-6d",abs_internal_energy);

        if ( model_is_ising )
            fprintf(file,"\nX = up-spin, 0 = down-spin");
        else  //  model is q-state Potts
            fprintf(file,"\nspin values: 1 (up) through %d",q);
        }
    }

fprintf(file,"\n");

if ( context == 2 )
    {
    if ( v_bonds )
        fprintf(file,"Virtual bonds created with probability = %.5f"
            " before cluster%s flipped.",v_bond_probability,
            ( dynamics_is_swendsen_wang ? "s" : "" ));
    else
        fprintf(file,"Actual bonds before cluster%s flipped.",
            ( dynamics_is_swendsen_wang ? "s" : "" ));
    }
else if ( context == 3 || context == 4 )
    {
    if ( v_bonds )
        fprintf(file,"Virtual bonds after cluster%s flipped.",
        ( dynamics_is_swendsen_wang ? "s" : "" ));
    else
        fprintf(file,"Actual bonds after cluster%s flipped.",
        ( dynamics_is_swendsen_wang ? "s" : "" ));

    if ( context == 4 )
        fprintf(file,"\nCluster grown from site i=%d,j=%d.  Cluster size = %d",
            w_i,w_j,ct_size_of_cluster);

    if ( dynamics_is_swendsen_wang || dynamics_is_wolff )
        {
        fprintf(file,"\nnum_clusters_traced_this_sweep = %u",num_clusters_traced_this_sweep);
        fprintf(file,"\nnum_spins_visited_this_sweep = %u",num_spins_visited_this_sweep);
        fprintf(file,"\nnum_spins_flipped_this_sweep = %u",num_spins_flipped_this_sweep);
        write_stack_data(file);
        }
    }

fprintf(file,"\n\n");

fprintf(file,"   ");
for ( i=0; i<size; i++ )
    fprintf(file,"%2d ",i);
fprintf(file,"\n");

for ( i=0; i<3*size; i++ )
    {
    if ( print_contracted_map )
        {
        //  Contract vertical double bonds to single bonds
        //  (if sure that the map is OK).  But if there are
        //  upper diagonal bonds then don't contract the map.
        if ( i>0 && !(i%3) )
            continue;
        }
    fprintf(file, "%2s ",( (i+2)%3 ? "  " : itoa(i/3,temp,10) ) );
    for ( j=0; j<3*size; j++ )
        fputc(slice[i][j],file);
    fprintf(file, " %-2s",( (i+2)%3 ? "  " : itoa(i/3,temp,10) ) );
    fputc('\n',file);
    }

fprintf(file,"   ");
for ( i=0; i<size; i++ )
    fprintf(file,"%2d ",i);
fprintf(file,"\n");

if ( model_is_q_state_potts && magn < 1.0 )
    {
    fprintf(file,"\nNumbers of q-state Potts spin values:");
    for ( j=0; j<2; j++ )
        {
        if ( 16*j + 1 > q )
            break;
        fprintf(file,"\n");
        for ( i=16*j; i<16*(j+1); i++ )
            {
            if ( i + 1 > q )
                break;
            fprintf(file,"%5d",i+1);
            }
        if ( q <= 16 || j == 1)
            fprintf(file," total");
        fprintf(file,"\n");
        for ( i=16*j; i<16*(j+1); i++ )
            {
            if ( i + 1 > q )
                break;
            fprintf(file,"%5d",q_state_potts_spins[i+1]);
            }
        }
    fprintf(file,"%5d",total_q_state_potts_spins);
    }

fprintf(file,"\n");
}

/*--------------------------------*/
static void write_3d_map(FILE *file,
                         int context,
                         int v_bonds)
{
}

/*--------------------------------*/
static void write_4d_map(FILE *file,
                         int context,
                         int v_bonds)
{
}