From a11bbf629d88379ecfdc68a6546bfc86b8403eaf Mon Sep 17 00:00:00 2001
From: Anthony Scemama <scemama@irsamc.ups-tlse.fr>
Date: Mon, 10 Jan 2022 13:19:22 +0100
Subject: [PATCH] Fix pseudos

---
 devel/trexio/export_trexio.irp.f | 128 ++++++++++++++++---------------
 1 file changed, 65 insertions(+), 63 deletions(-)

diff --git a/devel/trexio/export_trexio.irp.f b/devel/trexio/export_trexio.irp.f
index 51ec111..2e3e66a 100644
--- a/devel/trexio/export_trexio.irp.f
+++ b/devel/trexio/export_trexio.irp.f
@@ -71,99 +71,101 @@ subroutine run
 ! Pseudo-potentials
 ! -----------------
 
-  print *, 'ECP'
+  if (do_pseudo) then
 
-  integer                       :: num
+    print *, 'ECP'
+    integer                       :: num
 
-  num = 0
-  do k=1,pseudo_klocmax
-    do i=1,nucl_num
-      if (pseudo_dz_k(i,k) /= 0.d0) then
-        num = num+1
-      end if
-    end do
-  end do
-
-  do l=0,pseudo_lmax
-    do k=1,pseudo_kmax
+    num = 0
+    do k=1,pseudo_klocmax
       do i=1,nucl_num
-        if (pseudo_dz_kl(i,k,l) /= 0.d0) then
+        if (pseudo_dz_k(i,k) /= 0.d0) then
           num = num+1
         end if
       end do
     end do
-  end do
 
-  integer, allocatable          :: ang_mom(:), nucleus_index(:), power(:), lmax(:)
-  double precision, allocatable :: exponent(:), coefficient(:)
-
-  allocate(ang_mom(num), nucleus_index(num), exponent(num), coefficient(num), power(num), &
-           lmax(nucl_num) )
-
-  do i=1,nucl_num
-    lmax(i) = -1
     do l=0,pseudo_lmax
       do k=1,pseudo_kmax
-        if (pseudo_dz_kl_transp(k,l,i) /= 0.d0) then
-          lmax(i) = max(lmax(i), l)
-        end if
+        do i=1,nucl_num
+          if (pseudo_dz_kl(i,k,l) /= 0.d0) then
+            num = num+1
+          end if
+        end do
       end do
     end do
-  end do
 
-  j = 0
-  do i=1,nucl_num
-    do k=1,pseudo_klocmax
-      if (pseudo_dz_k_transp(k,i) /= 0.d0) then
-        j = j+1
-        ang_mom(j) = lmax(i)+1
-        nucleus_index(j) = i
-        exponent(j) = pseudo_dz_k_transp(k,i)
-        coefficient(j) = pseudo_v_k_transp(k,i)
-        power(j) = pseudo_n_k_transp(k,i)
-      end if
+    integer, allocatable          :: ang_mom(:), nucleus_index(:), power(:), lmax(:)
+    double precision, allocatable :: exponent(:), coefficient(:)
+
+    allocate(ang_mom(num), nucleus_index(num), exponent(num), coefficient(num), power(num), &
+            lmax(nucl_num) )
+
+    do i=1,nucl_num
+      lmax(i) = -1
+      do l=0,pseudo_lmax
+        do k=1,pseudo_kmax
+          if (pseudo_dz_kl_transp(k,l,i) /= 0.d0) then
+            lmax(i) = max(lmax(i), l)
+          end if
+        end do
+      end do
     end do
 
-    do l=0,lmax(i)
-      do k=1,pseudo_kmax
-        if (pseudo_dz_kl_transp(k,l,i) /= 0.d0) then
+    j = 0
+    do i=1,nucl_num
+      do k=1,pseudo_klocmax
+        if (pseudo_dz_k_transp(k,i) /= 0.d0) then
           j = j+1
-          ang_mom(j) = l
+          ang_mom(j) = lmax(i)+1
           nucleus_index(j) = i
-          exponent(j) = pseudo_dz_kl_transp(k,l,i)
-          coefficient(j) = pseudo_v_kl_transp(k,l,i)
-          power(j) = pseudo_n_kl_transp(k,l,i)
+          exponent(j) = pseudo_dz_k_transp(k,i)
+          coefficient(j) = pseudo_v_k_transp(k,i)
+          power(j) = pseudo_n_k_transp(k,i)
         end if
       end do
+
+      do l=0,lmax(i)
+        do k=1,pseudo_kmax
+          if (pseudo_dz_kl_transp(k,l,i) /= 0.d0) then
+            j = j+1
+            ang_mom(j) = l
+            nucleus_index(j) = i
+            exponent(j) = pseudo_dz_kl_transp(k,l,i)
+            coefficient(j) = pseudo_v_kl_transp(k,l,i)
+            power(j) = pseudo_n_kl_transp(k,l,i)
+          end if
+        end do
+      end do
     end do
-  end do
 
 
-  lmax(:) = lmax(:)+1
-  rc = trexio_write_ecp_max_ang_mom_plus_1(f, lmax)
-  call check_success(rc)
+    lmax(:) = lmax(:)+1
+    rc = trexio_write_ecp_max_ang_mom_plus_1(f, lmax)
+    call check_success(rc)
 
-  rc = trexio_write_ecp_z_core(f, int(nucl_charge_remove))
-  call check_success(rc)
+    rc = trexio_write_ecp_z_core(f, int(nucl_charge_remove))
+    call check_success(rc)
 
-  rc = trexio_write_ecp_num(f, num)
-  call check_success(rc)
+    rc = trexio_write_ecp_num(f, num)
+    call check_success(rc)
 
-  rc = trexio_write_ecp_ang_mom(f, ang_mom)
-  call check_success(rc)
+    rc = trexio_write_ecp_ang_mom(f, ang_mom)
+    call check_success(rc)
 
-  rc = trexio_write_ecp_nucleus_index(f, nucleus_index)
-  call check_success(rc)
+    rc = trexio_write_ecp_nucleus_index(f, nucleus_index)
+    call check_success(rc)
 
-  rc = trexio_write_ecp_exponent(f, exponent)
-  call check_success(rc)
+    rc = trexio_write_ecp_exponent(f, exponent)
+    call check_success(rc)
 
-  rc = trexio_write_ecp_coefficient(f, coefficient)
-  call check_success(rc)
+    rc = trexio_write_ecp_coefficient(f, coefficient)
+    call check_success(rc)
 
-  rc = trexio_write_ecp_power(f, power)
-  call check_success(rc)
+    rc = trexio_write_ecp_power(f, power)
+    call check_success(rc)
 
+  endif
 
 
 ! Basis