00001 subroutine invhij
00002 use grid_parameter, only : nrg, ntg, npg
00003 use def_metric_hij, only : hxxd, hxyd, hxzd, hyyd, hyzd, hzzd, &
00004 & hxxu, hxyu, hxzu, hyyu, hyzu, hzzu
00005 implicit none
00006 real(8) :: hxx, hxy, hxz, hyx, hyy, hyz, hzx, hzy, hzz,
00007 hod1, hod2, hod3, detgm, detgmi, &
00008 gmxxu, gmxyu, gmxzu, gmyxu, gmyyu, gmyzu, &
00009 gmzxu, gmzyu, gmzzu
00010 integer :: ipg, itg, irg
00011
00012 do ipg = 0, npg
00013 do itg = 0, ntg
00014 do irg = 0, nrg
00015
00016 hxx = hxxd(irg,itg,ipg)
00017 hxy = hxyd(irg,itg,ipg)
00018 hxz = hxzd(irg,itg,ipg)
00019 hyx = hxy
00020 hyy = hyyd(irg,itg,ipg)
00021 hyz = hyzd(irg,itg,ipg)
00022 hzx = hxz
00023 hzy = hyz
00024 hzz = hzzd(irg,itg,ipg)
00025
00026 hod1 = hxx + hyy + hzz
00027 hod2 = hxx*hyy + hxx*hzz + hyy*hzz &
00028 & - hxy*hyx - hxz*hzx - hyz*hzy
00029 hod3 = hxx*hyy*hzz + hxy*hyz*hzx + hxz*hyx*hzy &
00030 & - hxx*hyz*hzy - hxy*hyx*hzz - hxz*hyy*hzx
00031 detgm = 1.0d0 + hod1 + hod2 + hod3
00032 detgmi = 1.0d0/detgm
00033
00034 hod1 = + hyy + hzz
00035 hod2 = + hyy*hzz - hyz*hzy
00036 gmxxu = (1.0d0 + hod1 + hod2)*detgmi
00037 hod1 = - hxy
00038 hod2 = + hxz*hzy - hxy*hzz
00039 gmxyu = (hod1 + hod2)*detgmi
00040 hod1 = - hxz
00041 hod2 = + hxy*hyz - hxz*hyy
00042 gmxzu = (hod1 + hod2)*detgmi
00043 hod1 = + hxx + hzz
00044 hod2 = + hxx*hzz - hxz*hzx
00045 gmyyu = (1.0d0 + hod1 + hod2)*detgmi
00046 hod1 = - hyz
00047 hod2 = + hxz*hyx - hxx*hyz
00048 gmyzu = (hod1 + hod2)*detgmi
00049 hod1 = + hxx + hyy
00050 hod2 = + hxx*hyy - hxy*hyx
00051 gmzzu = (1.0d0 + hod1 + hod2)*detgmi
00052 gmyxu = gmxyu
00053 gmzxu = gmxzu
00054 gmzyu = gmyzu
00055
00056 hxxu(irg,itg,ipg) = gmxxu - 1.0d0
00057 hxyu(irg,itg,ipg) = gmxyu
00058 hxzu(irg,itg,ipg) = gmxzu
00059 hyyu(irg,itg,ipg) = gmyyu - 1.0d0
00060 hyzu(irg,itg,ipg) = gmyzu
00061 hzzu(irg,itg,ipg) = gmzzu - 1.0d0
00062
00063 end do
00064 end do
00065 end do
00066
00067 end subroutine invhij
