00001 subroutine liegmab_gridpoint
00002   use grid_parameter, only : nrg, ntg, npg
00003   use def_Lie_derivatives_grid, only : elpxx_grid, elpxy_grid, elpxz_grid, &
00004   &                               elpyy_grid, elpyz_grid, elpzz_grid, &
00005   &                               rlpxx_grid, rlpxy_grid, rlpxz_grid, &
00006   &                               rlpyy_grid, rlpyz_grid, rlpzz_grid
00007   use def_metric_hij, only : hxxd, hxyd, hxzd, hyyd, hyzd, hzzd, &
00008   &                          hxxu, hxyu, hxzu, hyyu, hyzu, hzzu
00009   use def_dvphi, only : dphiu
00010   use interface_grdphi_gridpoint_type0
00011   implicit none
00012   real(8) :: dhxxdp, dhxxdx, dhxxdy, dhxxdz, 
00013             dhxydp, dhxydx, dhxydy, dhxydz, &
00014             dhxzdp, dhxzdx, dhxzdy, dhxzdz, &
00015             dhyxdp, dhyxdx, dhyxdy, dhyxdz, &
00016             dhyydp, dhyydx, dhyydy, dhyydz, &
00017             dhyzdp, dhyzdx, dhyzdy, dhyzdz, &
00018             dhzxdp, dhzxdx, dhzxdy, dhzxdz, &
00019             dhzydp, dhzydx, dhzydy, dhzydz, &
00020             dhzzdp, dhzzdx, dhzzdy, dhzzdz, &
00021             fa13, trlie, &
00022             gmxxd, gmxyd, gmxzd, gmxxu, gmxyu, gmxzu, &
00023             gmyxd, gmyyd, gmyzd, gmyxu, gmyyu, gmyzu, &
00024             gmzxd, gmzyd, gmzzd, gmzxu, gmzyu, gmzzu
00025   integer :: ipg, irg, itg
00026 
00027 
00028 
00029   fa13 = 1.0d0/3.0d0
00030   dphiu(1:3,1:3) = 0.0d0
00031   dphiu(1,2) =-1.0d0
00032   dphiu(2,1) = 1.0d0
00033 
00034   do ipg = 0, npg
00035     do itg = 0, ntg
00036       do irg = 0, nrg
00037 
00038         gmxxd=hxxd(irg,itg,ipg)
00039         gmxyd=hxyd(irg,itg,ipg)
00040         gmxzd=hxzd(irg,itg,ipg)
00041         gmyyd=hyyd(irg,itg,ipg)
00042         gmyzd=hyzd(irg,itg,ipg)
00043         gmzzd=hzzd(irg,itg,ipg)
00044         gmxxd = gmxxd + 1.0d0
00045         gmyyd = gmyyd + 1.0d0
00046         gmzzd = gmzzd + 1.0d0
00047         gmyxd = gmxyd
00048         gmzxd = gmxzd
00049         gmzyd = gmyzd
00050         gmxxu=hxxu(irg,itg,ipg)
00051         gmxyu=hxyu(irg,itg,ipg)
00052         gmxzu=hxzu(irg,itg,ipg)
00053         gmyyu=hyyu(irg,itg,ipg)
00054         gmyzu=hyzu(irg,itg,ipg)
00055         gmzzu=hzzu(irg,itg,ipg)
00056         gmxxu = gmxxu + 1.0d0
00057         gmyyu = gmyyu + 1.0d0
00058         gmzzu = gmzzu + 1.0d0
00059         gmyxu = gmxyu
00060         gmzxu = gmxzu
00061         gmzyu = gmyzu
00062 
00063         call grdphi_gridpoint_type0(hxxd,dhxxdp,irg,itg,ipg)
00064         call grdphi_gridpoint_type0(hxyd,dhxydp,irg,itg,ipg)
00065         call grdphi_gridpoint_type0(hxzd,dhxzdp,irg,itg,ipg)
00066         call grdphi_gridpoint_type0(hyyd,dhyydp,irg,itg,ipg)
00067         call grdphi_gridpoint_type0(hyzd,dhyzdp,irg,itg,ipg)
00068         call grdphi_gridpoint_type0(hzzd,dhzzdp,irg,itg,ipg)
00069 
00070         rlpxx_grid(irg,itg,ipg) = dhxxdp + gmxyd*dphiu(2,1) + gmyxd*dphiu(2,1)
00071         rlpxy_grid(irg,itg,ipg) = dhxydp + gmxxd*dphiu(1,2) + gmyyd*dphiu(2,1)
00072         rlpxz_grid(irg,itg,ipg) = dhxzdp
00073         rlpyy_grid(irg,itg,ipg) = dhyydp + gmyxd*dphiu(1,2) + gmxyd*dphiu(1,2)
00074         rlpyz_grid(irg,itg,ipg) = dhyzdp
00075         rlpzz_grid(irg,itg,ipg) = dhzzdp
00076 
00077         dhxxdp = rlpxx_grid(irg,itg,ipg)
00078         dhxydp = rlpxy_grid(irg,itg,ipg)
00079         dhxzdp = rlpxz_grid(irg,itg,ipg)
00080         dhyydp = rlpyy_grid(irg,itg,ipg)
00081         dhyzdp = rlpyz_grid(irg,itg,ipg)
00082         dhzzdp = rlpzz_grid(irg,itg,ipg)
00083         dhyxdp = dhxydp
00084         dhzxdp = dhxzdp
00085         dhzydp = dhyzdp
00086 
00087         trlie = gmxxu*dhxxdp + gmxyu*dhxydp + gmxzu*dhxzdp &
00088            &      + gmyxu*dhyxdp + gmyyu*dhyydp + gmyzu*dhyzdp &
00089            &      + gmzxu*dhzxdp + gmzyu*dhzydp + gmzzu*dhzzdp
00090 
00091         elpxx_grid(irg,itg,ipg) = dhxxdp - fa13*gmxxd*trlie
00092         elpxy_grid(irg,itg,ipg) = dhxydp - fa13*gmxyd*trlie
00093         elpxz_grid(irg,itg,ipg) = dhxzdp - fa13*gmxzd*trlie
00094         elpyy_grid(irg,itg,ipg) = dhyydp - fa13*gmyyd*trlie
00095         elpyz_grid(irg,itg,ipg) = dhyzdp - fa13*gmyzd*trlie
00096         elpzz_grid(irg,itg,ipg) = dhzzdp - fa13*gmzzd*trlie
00097 
00098       end do
00099     end do
00100   end do
00101 
00102 end subroutine liegmab_gridpoint