fullpos3D.xmds

<?xmds version="1.0"?>
<!--Stochastic Superchemistry simulation-->

<simulation>

  <stochastic>yes</stochastic> 
  <prop_dim>t</prop_dim>
  <error_check>yes</error_check>
  <paths>1</paths>
  <use_mpi>no</use_mpi> 
  <seed>37 41</seed>
  <noises>8</noises>

  <globals>
    <![CDATA[
    const double noise = 1.0;
    const double hbar = 1.05500000000e-34;
    const double M = 1.409539200000000e-25;
    const double omegax = 0.58976353090742;
    const double omegay = 0.58976353090742;
    const double omegaz = 0.58976353090742/30;
    const double U11 = 2.974797272874263e-51;
    const double U13 = -1.417820412490823e-50;
    const double U33 = 2.974797272874263e-51;
    const double inum = 1.0e6;
    const double Uoh11 = U11/hbar;
    const double Uoh13 = U13/hbar;
    const double Uoh33 = U33/hbar;
    const double mu = pow(15*inum*U11*omegax*omegay*omegaz/M_PI/4,0.4)*pow(M,0.6)/2;
    const double delta = 1.0e9;
    const double F = 2.0e-2;
    const double g = sqrt(Uoh11*2/delta);
    const double loss11=1.0e-2;
    const double loss12=1.6e-22;
    const double loss31=1.0e-2;
    const double loss32=1.6e-22;
    const double loss132=8.0e-17;
    const double chi = F*g*delta;
    const double biggamma = g*g*delta/2;
    const double gam13 = Uoh13/chi;
    const double gam33 = Uoh33/chi;
    const double gameff = (Uoh11-biggamma)/chi;
    const double gamloss11=loss11/2/chi;
    const double gamloss12=loss12/chi;
    const double gamloss31=loss31/2/chi;
    const double gamloss32=loss32/chi;
    const double gamloss132=loss132/chi;
    const double cnoise =  noise/sqrt(2); 
    ]]>
  </globals>  

  <field>
    <dimensions>x y z</dimensions>
    <lattice>32 32 16</lattice>
    <domains>(-1.2e-4,1.2e-4) (-1.2e-4,1.2e-4) (-8.0e-3,8.0e-3)</domains>
    <samples> 1 1 </samples>

    <vector>
      <name> vc1 </name>
      <type>double</type>
      <components>vcore V1r V3r gV1r gV3r</components>
      <fourier_space>no no no</fourier_space>
      <![CDATA[
      vcore = (omegax*omegax*x*x+omegay*omegay*y*y+omegaz*omegaz*z*z);
      V1r = 0.5*M*vcore/hbar/chi -(gameff+gam13/2)/2/(dx*dy*dz);
      V3r = M*vcore/hbar/chi -(gam13/2+gam33)/2/(dx*dy*dz);
      gV1r = 0.5*M*vcore/hbar/chi;
      gV3r = M*vcore/hbar/chi;
      ]]>
    </vector>

    <vector>
      <name> main </name>
      <type>complex</type>
      <components>phi1a phi1b phi3a phi3b gphi1a gphi3a</components>
      <fourier_space>no no no</fourier_space>
      <vectors> vc1 </vectors>
      <![CDATA[
      const double realfn = (mu-0.5*M*vcore)/Uoh11/hbar;

      phi1a = realfn>0. ? complex(sqrt(realfn),0) : complex(0,0);
      phi1b = realfn>0. ? complex(sqrt(realfn),0) : complex(0,0);
      phi3a = complex(0,0);
      phi3b = complex(0,0);
      gphi1a = realfn>0. ? complex(sqrt(realfn),0) : complex(0,0);
      gphi3a = complex(0,0);
      ]]>
    </vector>
  </field>

  <sequence>
    <integrate>
      <algorithm>RK4IP</algorithm>
      <interval>1e-10</interval>
      <lattice>5</lattice>
      <samples>1 5</samples>
      <k_operators>
	<constant>yes</constant>
	<operator_names> L2p L2n L4p L4n </operator_names>
	<![CDATA[
	L2p = complex(0,-hbar/M/2/chi*(kx*kx+ky*ky+kz*kz));
	L2n = complex(0, hbar/M/2/chi*(kx*kx+ky*ky+kz*kz));
	L4p = complex(0,-hbar/M/4/chi*(kx*kx+ky*ky+kz*kz));
	L4n = complex(0, hbar/M/4/chi*(kx*kx+ky*ky+kz*kz));
	]]>
      </k_operators>
      <vectors> main vc1 </vectors>
      <![CDATA[
      const complex dens1 = phi1b*phi1a;
      const complex dens3 = phi3b*phi3a;

      const double gdens1 = (gphi1a.re*gphi1a.re+gphi1a.im*gphi1a.im);
      const double gdens3 = (gphi3a.re*gphi3a.re+gphi3a.im*gphi3a.im);

      dphi1a_dt = L2p[phi1a] + (-i*V1r-gamloss11+(gamloss132/2+gamloss12)/2/(dx*dy*dz))*phi1a 
      +(-i*gameff-gamloss12)*dens1*phi1a 
      -(i*gam13+gamloss132)*dens3*phi1a -i*phi1b*phi3a 
      +c_sqrt(-i*phi3a-(i*gameff+gamloss12)*phi1a*phi1a)*noise*n_1
      -c_sqrt((i*gam13+gamloss132)*phi1a*phi3a)*cnoise*(n_5-i*n_7);

      dphi1b_dt = L2n[phi1b] + (i*V1r-gamloss11+(gamloss132/2+gamloss12)/2/(dx*dy*dz))*phi1b 
      +(i*gameff-gamloss12)*dens1*phi1b 
      +(i*gam13-gamloss132)*dens3*phi1b +i*phi1a*phi3b
      +c_sqrt(i*phi3b+(i*gameff-gamloss12)*phi1b*phi1b)*noise*n_2
      +c_sqrt((gamloss132-i*gam13)*phi1b*phi3b)*cnoise*(n_6+i*n_8);

      dphi3a_dt = L4p[phi3a] + (-i*V3r-gamloss31+(gamloss132/2+gamloss32)/2/(dx*dy*dz))*phi3a 
      +(-i*gam33-gamloss32)*dens3*phi3a 
      -i*0.5*phi1a*phi1a -(i*gam13+gamloss132)*dens1*phi3a
      +c_sqrt(-i*gam33-gamloss32)*phi3a*noise*n_3
      +c_sqrt((gamloss132+i*gam13)*phi1a*phi3a)*cnoise*(n_5+i*n_7);

      dphi3b_dt = L4n[phi3b] + (i*V3r-gamloss31+(gamloss132/2+gamloss32)/2/(dx*dy*dz))*phi3b 
      +(i*gam33-gamloss32)*dens3*phi3b 
      + i*0.5*phi1b*phi1b +(i*gam13-gamloss132)*dens1*phi3b
      +c_sqrt(i*gam33-gamloss32)*phi3b*noise*n_4
      +c_sqrt((gamloss132-i*gam13)*phi1b*phi3b)*cnoise*(n_6-i*n_8);
      
      dgphi1a_dt = L2p[gphi1a] + (-i*gV1r-gamloss11)*gphi1a +(-i*gameff-gamloss12)*gdens1*gphi1a 
      -(i*gam13+gamloss132)*gdens3*gphi1a-i*conj(gphi1a)*gphi3a;

      dgphi3a_dt = L4p[gphi3a] + (-i*gV3r-gamloss31)*gphi3a +(-i*gam33-gamloss32)*gdens3*gphi3a 
      -i*0.5*gphi1a*gphi1a +(i*gam13-gamloss132)*gdens1*gphi3a;
      ]]>
    </integrate>
  </sequence>

  <output>
    <group>
      <sampling>
	<fourier_space>   no     no no</fourier_space>
	<lattice>         16    4  4</lattice>
	<moments>atoms molecules gatoms gmolecules</moments>
	<![CDATA[
	atoms=phi1b*phi1a;
	molecules=phi3b*phi3a;
	gatoms=conj(gphi1a)*gphi1a;
	gmolecules=conj(gphi3a)*gphi3a;
	]]>
      </sampling>
    </group>

    <group>
      <sampling>
	<fourier_space>    no   no  no</fourier_space>
	<lattice>           0    0   0</lattice>
	<moments>rn_1 rn_2 grn_1 grn_2 excitedn</moments>
	<![CDATA[
	rn_1 = phi1b*phi1a;
	rn_2 = phi3b*phi3a;
	grn_1 = conj(gphi1a)*gphi1a;
	grn_2 = conj(gphi3a)*gphi3a;
	excitedn = g*g/4*phi1b*phi1b*phi1a*phi1a+F*F*phi3b*phi3a
	- F*g/2*(phi1b*phi1b*phi3a+phi1a*phi1a*phi3b);
	]]>
      </sampling>
    </group>
  </output>
</simulation>