FirstMultiply.cpp File Reference

Multiplication \( v_0 = H v_1 \) at the first step for TPQ mode ( \( v_1 \) is the random or inputted vector). More...

#include "FirstMultiply.hpp"
#include "expec_energy_flct.hpp"
#include "common/setmemory.hpp"
#include "wrapperMPI.hpp"
#include "CalcTime.hpp"
#include "mltplyCommon.hpp"
#include "expec_cisajs.hpp"
#include "expec_cisajscktaltdc.hpp"

Go to the source code of this file.

Functions
int	FirstMultiply (struct BindStruct *X)
	Multiplication \( v_0 = H v_1 \) at the first step for TPQ mode ( \( v_1 \) is the random or inputted vector). More...

Detailed Description

Multiplication \( v_0 = H v_1 \) at the first step for TPQ mode ( \( v_1 \) is the random or inputted vector).

Author

Takahiro Misawa (The University of Tokyo)

Kazuyoshi Yoshimi (The University of Tokyo)

Version

0.1

Definition in file FirstMultiply.cpp.

Function Documentation

FirstMultiply()

int FirstMultiply

(

struct BindStruct *

X

)

Multiplication \( v_0 = H v_1 \) at the first step for TPQ mode ( \( v_1 \) is the random or inputted vector).

Parameters

rand_i	[in] A rundom number seed for giving the initial vector \( v_1 \).
X	[in] Struct to get information of the vector \( v_1 \) for the first step calculation.

Return values

-1	fail the multiplication \( v_0 = H v_1 \).
0	succeed the multiplication \( v_0 = H v_1 \).

Version

0.1

Author

Takahiro Misawa (The University of Tokyo)

Kazuyoshi Yoshimi (The University of Tokyo)

Compute expectation value at infinite temperature

Definition at line 42 of file FirstMultiply.cpp.

Referenced by CalcByTPQ().

                                        {

  long int i, i_max;
  double dnorm;
  double Ns;
  long int u_long_i;
  dsfmt_t dsfmt;
  int mythread, rand_i, iret;

  Ns = 1.0*X->Def.NsiteMPI;
  i_max = X->Check.idim_max;

  for (rand_i = 0; rand_i < NumAve; rand_i++) {
#pragma omp parallel default(none) private(i, mythread, u_long_i, dsfmt) \
shared(I, v0, v1, nthreads, myrank, rand_i, X, stdoutMPI) \
firstprivate(i_max)
  {
#pragma omp for
    for (i = 1; i <= i_max; i++) {
      v0[i][rand_i] = 0.0;
      v1[i][rand_i] = 0.0;
    }
    /*
    Initialise MT
    */
#ifdef _OPENMP
    mythread = omp_get_thread_num();
#else
    mythread = 0;
#endif
    u_long_i = 123432 + (rand_i + 1)*labs(X->Def.initial_iv) + mythread + nthreads * myrank;
    dsfmt_init_gen_rand(&dsfmt, u_long_i);

    if (X->Def.iInitialVecType == 0) {

      StartTimer(3101);
#pragma omp for
      for (i = 1; i <= i_max; i++)
        v1[i][rand_i] = 2.0*(dsfmt_genrand_close_open(&dsfmt) - 0.5) + 2.0*(dsfmt_genrand_close_open(&dsfmt) - 0.5)*I;
      }/*if (X->Def.iInitialVecType == 0)*/
      else {
#pragma omp for
        for (i = 1; i <= i_max; i++)
          v1[i][rand_i] = 2.0*(dsfmt_genrand_close_open(&dsfmt) - 0.5);
      }
      StopTimer(3101);
    }/*#pragma omp parallel*/
    /*
    Normalize v
    */
  dnorm = 0.0;
#pragma omp parallel for default(none) private(i) \
shared(v1, i_max, rand_i) reduction(+:dnorm)
      for (i = 1; i <= i_max; i++)
        dnorm += real(conj(v1[i][rand_i])*v1[i][rand_i]);
    dnorm = SumMPI_d(dnorm);
    dnorm = sqrt(dnorm);
    global_1st_norm[rand_i] = dnorm;
#pragma omp parallel for default(none) private(i) shared(v1,rand_i) firstprivate(i_max, dnorm)
    for (i = 1; i <= i_max; i++) v1[i][rand_i] = v1[i][rand_i] / dnorm;
  }/*for (rand_i = 0; rand_i < NumAve; rand_i++)*/

  TimeKeeperWithRandAndStep(X, "%s_TimeKeeper.dat", "set %d step %d:TPQ begins: %s", "a", rand_i, step_i);
  X->Def.istep = 0;
  StartTimer(3300);
  iret = expec_cisajs(X, NumAve, v0, v1);
  StopTimer(3300);
  if (iret != 0) return -1;

  StartTimer(3400);
  iret = expec_cisajscktaltdc(X, NumAve, v0, v1);
  StopTimer(3400);
  if (iret != 0) return -1;

#pragma omp parallel for default(none) private(i,rand_i) shared(v0,v1,i_max,NumAve)
  for (i = 1; i <= i_max; i++) 
    for (rand_i = 0; rand_i < NumAve; rand_i++) v0[i][rand_i] = v1[i][rand_i];
  StartTimer(3102);
  if(expec_energy_flct(X, NumAve, v0, v1) !=0){
    StopTimer(3102);
    return -1;
  }
  StopTimer(3102);

  for (rand_i = 0; rand_i < NumAve; rand_i++) {
#pragma omp parallel for default(none) private(i) shared(v0, v1, list_1,rand_i) \
firstprivate(i_max, Ns, LargeValue, myrank)
    for (i = 1; i <= i_max; i++) {
      v0[i][rand_i] = LargeValue * v1[i][rand_i] - v0[i][rand_i] / Ns;
    }

    dnorm = 0.0;
#pragma omp parallel for default(none) private(i) shared(v0,rand_i) \
firstprivate(i_max) reduction(+:dnorm)
    for (i = 1; i <= i_max; i++)
      dnorm += real(conj(v0[i][rand_i])*v0[i][rand_i]);
    dnorm = SumMPI_d(dnorm);
    dnorm = sqrt(dnorm);
    global_norm[rand_i] = dnorm;
#pragma omp parallel for default(none) private(i) shared(v0,rand_i) firstprivate(i_max, dnorm)
    for (i = 1; i <= i_max; i++) v0[i][rand_i] = v0[i][rand_i] / dnorm;
  }/*for (rand_i = 0; rand_i < NumAve; rand_i++)*/
  TimeKeeperWithRandAndStep(X, "%s_TimeKeeper.dat", "set %d step %d:TPQ finishes: %s", "a", rand_i, step_i);
  return 0;
}