Not working

src/AO/hessian.irp.f

@@ -0,0 +1,119 @@
+BEGIN_PROVIDER [ real, ao_axis_hess, (ao_num,elec_num,3,3) ] 
+ implicit none 
+ include 'types.F' 
+ 
+ BEGIN_DOC 
+! Gradients of the cartesian polynomial part of the atomic orbitals. 
+ END_DOC 
+ integer :: i, j, l 
+ real,save   :: real_of_int(-2:10) 
+ data real_of_int /0.,0.,0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10./ 
+ 
+ do j=1,elec_num 
+  do i=1,ao_num 
+!       a*b*(x-x0)**(a-1)*(y-y0)**(b-1)*(z-z0)**c 
+    ao_axis_hess(i,j,1,1) = real_of_int(ao_power(i,1)) * real_of_int(ao_power(i,1)-1) & 
+                              * ao_axis_power( ao_power(i,1)-2, 1 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,2)  , 2 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,3)  , 3 , ao_nucl(i) , j ) 
+ 
+    ao_axis_hess(i,j,2,1) = real_of_int(ao_power(i,1)) * real_of_int(ao_power(i,2)) & 
+                              * ao_axis_power( ao_power(i,1)-1, 1 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,2)-1, 2 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,3)  , 3 , ao_nucl(i) , j ) 
+ 
+    ao_axis_hess(i,j,3,1) = real_of_int(ao_power(i,1)) * real_of_int(ao_power(i,3)) & 
+                              * ao_axis_power( ao_power(i,1)-1, 1 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,2)  , 2 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,3)-1, 3 , ao_nucl(i) , j ) 
+ 
+    ao_axis_hess(i,j,1,2) = ao_axis_hess(i,j,2,1) 
+    ao_axis_hess(i,j,2,2) = real_of_int(ao_power(i,2)) * real_of_int(ao_power(i,2)-1) & 
+                              * ao_axis_power( ao_power(i,1)  , 1 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,2)-2, 2 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,3)  , 3 , ao_nucl(i) , j ) 
+ 
+    ao_axis_hess(i,j,3,2) = real_of_int(ao_power(i,2)) * real_of_int(ao_power(i,3)) & 
+                              * ao_axis_power( ao_power(i,1)  , 1 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,2)-1, 2 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,3)-1, 3 , ao_nucl(i) , j ) 
+ 
+    ao_axis_hess(i,j,1,3) = ao_axis_hess(i,j,3,1) 
+    ao_axis_hess(i,j,2,3) = ao_axis_hess(i,j,3,2) 
+    ao_axis_hess(i,j,3,3) = real_of_int(ao_power(i,3)) * real_of_int(ao_power(i,3)-1) & 
+                              * ao_axis_power( ao_power(i,1)  , 1 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,2)  , 2 , ao_nucl(i) , j ) & 
+                              * ao_axis_power( ao_power(i,3)-2, 3 , ao_nucl(i) , j ) 
+  enddo 
+ enddo 
+ 
+END_PROVIDER 
+
+BEGIN_PROVIDER [ real, ao_axis_lapl_grad, (ao_num,elec_num,3) ] 
+ implicit none 
+ include 'types.F' 
+ BEGIN_DOC 
+! Laplacian of the cartesian atomic orbitals  
+ END_DOC 
+ integer :: i, j, l 
+ 
+ do j=1,elec_num 
+ 
+   do i=1,ao_num 
+     real, save :: real_of_int(-3:10) 
+     data real_of_int /0.,0.,0.,0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10./ 
+ 
+      ao_axis_lapl_grad(i,j,1) = real_of_int(ao_power(i,1)-2) *  & 
+         real_of_int(ao_power(i,1)-1) * real_of_int(ao_power(i,1)) * & 
+         ao_axis_power( ao_power(i,1)-3, 1 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,2), 2 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,3), 3 , ao_nucl(i) , j ) & 
+        + real_of_int(ao_power(i,1)) *  & 
+         real_of_int(ao_power(i,2)-1) * real_of_int(ao_power(i,2)) * & 
+         ao_axis_power( ao_power(i,1)-1, 1 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,2)-2, 2 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,3), 3 , ao_nucl(i) , j ) & 
+        + real_of_int(ao_power(i,1)) *  & 
+         real_of_int(ao_power(i,3)-1) * real_of_int(ao_power(i,3)) * & 
+         ao_axis_power( ao_power(i,1)-1, 1 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,2), 2 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,3)-2, 3 , ao_nucl(i) , j ) 
+ 
+      ao_axis_lapl_grad(i,j,2) = real_of_int(ao_power(i,2)-2) *  & 
+         real_of_int(ao_power(i,2)-1) * real_of_int(ao_power(i,2)) * & 
+         ao_axis_power( ao_power(i,2)-3, 2 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,3), 3 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,1), 1 , ao_nucl(i) , j ) & 
+        + real_of_int(ao_power(i,2)) *  & 
+         real_of_int(ao_power(i,3)-1) * real_of_int(ao_power(i,3)) * & 
+         ao_axis_power( ao_power(i,2)-1, 2 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,3)-2, 3 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,1), 1 , ao_nucl(i) , j ) & 
+        + real_of_int(ao_power(i,1)) *  & 
+         real_of_int(ao_power(i,1)-1) * real_of_int(ao_power(i,2)) * & 
+         ao_axis_power( ao_power(i,2)-1, 2 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,3), 3 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,1)-2, 1 , ao_nucl(i) , j ) 
+ 
+      ao_axis_lapl_grad(i,j,3) = real_of_int(ao_power(i,3)-2) *  & 
+         real_of_int(ao_power(i,3)-1) * real_of_int(ao_power(i,3)) * & 
+         ao_axis_power( ao_power(i,3)-3, 3 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,1), 1 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,2), 2 , ao_nucl(i) , j ) & 
+        + real_of_int(ao_power(i,3)) *  & 
+         real_of_int(ao_power(i,1)-1) * real_of_int(ao_power(i,1)) * & 
+         ao_axis_power( ao_power(i,3)-1, 3 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,1)-2, 1 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,2), 2 , ao_nucl(i) , j ) & 
+        + real_of_int(ao_power(i,2)) *  & 
+         real_of_int(ao_power(i,2)-1) * real_of_int(ao_power(i,3)) * & 
+         ao_axis_power( ao_power(i,3)-1, 3 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,1), 1 , ao_nucl(i) , j ) * & 
+         ao_axis_power( ao_power(i,2)-2, 2 , ao_nucl(i) , j ) 
+ 
+   enddo 
+ 
+ enddo 
+ 
+END_PROVIDER 
+