anderkve / FYS3150

https://anderkve.github.io/FYS3150
26 stars 14 forks source link

Phase transition #90

Closed JohanCarlsen closed 1 year ago

JohanCarlsen commented 1 year ago

Hi! So I've produced plots of $<\epsilon>$, $<|m|>$, $C_V$ and $\chi$ for different temperatures. I've used t=linspace(2.1, 2.4, 6) and run 1.000.000 MC cycles on the systems. From the compendium, I notice the shape of my curves are "swapped" for $<\epsilon>$ and $<|m|>$. Also, the plot for $<|m|>$ does not show the clear step-function shape I expected. Does anyone have an idea about what the reason for this might be? phase_transition

int numTempElements = 6;
    vec tempList = linspace(2.1, 2.4, numTempElements);

    mat phaseRes = mat(numTempElements, 4);

    #ifdef _OPENMP
    {
      #pragma omp parallel
      {
        const int threadID = omp_get_thread_num();

        mt19937 generator;
        unsigned int threadSeed = baseSeed + threadID;
        generator.seed(threadSeed);

        // Timing algorithm
        auto t1 = chrono::high_resolution_clock::now();

        #pragma omp for
        for (int t = 0; t < numTempElements; t++)
        {
          int burnIn = nCyclesPerThread / 10;

          double T_phase = tempList(t);
          double threadE = 0, threadM = 0;
          double sumE = 0, sumEE = 0, sumM = 0, sumMM = 0;
          double avgE = 0, avgEE = 0, avgM = 0, avgMM = 0;
          double avg_e = 0, avg_m = 0;
          double heatCap = 0, X = 0;

          mat threadLattice = initialize_lattice(L, T_phase, threadE, threadM, false);

          for (int burn = 0; burn < burnIn; burn++)
          {
            mcmc(threadLattice, generator, L, threadE, threadM, T_phase);
          }

          for (int n = burnIn; n < nCyclesPerThread; n++)
          {
            mcmc(threadLattice, generator, L, threadE, threadM, T_phase);

            sumE += threadE;
            sumEE += threadE * threadE;
            sumM += fabs(threadM);
            sumMM += threadM * threadM;
          }

          avg_e = sumE / (nCyclesPerThread - burnIn * N);
          avg_m = sumM / (nCyclesPerThread - burnIn * N);

          avgE = sumE / (nCyclesPerThread - burnIn);
          avgEE = sumEE / (nCyclesPerThread - burnIn);
          avgM = sumM / (nCyclesPerThread - burnIn);
          avgMM = sumMM / (nCyclesPerThread - burnIn);

          heatCap = (avgEE - avgE * avgE) / (N * T_phase * T_phase);
          X = (avgMM - avgM * avgM) / (N * T_phase);

          phaseRes(t, 0) = avg_e;
          phaseRes(t, 1) = avg_m;
          phaseRes(t, 2) = heatCap;
          phaseRes(t, 3) = X;
anderkve commented 1 year ago

Hi @JohanCarlsen!

I haven't looked closely at it, but one thing you can look into:

The values you get for <e> and <|m|> are strange. They are too large in magnitude (think about what the largest possible values should be). And how can you get a negative number for the expectation value <|m|> for the absolute magnetization? So I think a place to start is to double-check the parentheses in this part :)

          avg_e = sumE / (nCyclesPerThread - burnIn * N);
          avg_m = sumM / (nCyclesPerThread - burnIn * N);

Other than that, I'd suggest you search for a possible sign error in your code. Perhaps it's somewhere in the mcmc function, so that the returned threadE and threadM are wrong? (I'm just guessing here.)

JohanCarlsen commented 1 year ago

So I finally fixed the code. I have zero clue about what was making the errors, but it is fixed now! wide_phase_transition

anderkve commented 1 year ago

Hooray! :)