diff --git a/src/pseudo.irp.f b/src/pseudo.irp.f
index dd01e01..156c25a 100644
--- a/src/pseudo.irp.f
+++ b/src/pseudo.irp.f
@@ -333,31 +333,54 @@ real function ylm(l,m,x,y,z,r_inv)
       continue
 
     case(1)
+      ylm = ylm * sq3 * r_inv
       select case (m)
         case (-1)
-          ylm = ylm * sq3 * y * r_inv
+          ylm = ylm * y 
         case (0)                      
-          ylm = ylm * sq3 * z * r_inv
+          ylm = ylm * z
         case (1)                      
-          ylm = ylm * sq3 * x * r_inv
+          ylm = ylm * x
+      end select
+
+    case(2)
+      ylm = ylm * r_inv * r_inv * sqrt(5.)
+      select case (m)
+        case(-2) 
+          ylm =  ylm * sqrt(3.) * x * y 
+        case(-1) 
+          ylm =  ylm * sqrt(3.) * y * z 
+        case(0) 
+          ylm =  ylm * 0.5 * (2.*z*z - x*x - y*y)
+        case(1) 
+          ylm =  ylm * sqrt(3.) * z * x
+        case(2) 
+          ylm =  ylm * 0.5 * sqrt(3.) * (x*x - y*y) 
+      end select
+
+    case(3)
+      ylm = ylm * r_inv * r_inv * r_inv * sqrt(7.)
+      select case (m)
+        case(-3) 
+          ylm =  ylm * 0.25 * sqrt(10.) * (3.*x*x - y*y)
+        case(-2) 
+          ylm =  ylm * sqrt(15.) * x*y*z
+        case(-1) 
+          ylm =  ylm * 0.25*sqrt(6.) * y * (4.*z*z - x*x -y*y)
+        case(0) 
+          ylm =  ylm * 0.5 * z * (2.*z*z - 3.*(x*x + y*y))
+        case(1) 
+          ylm =  ylm * 0.25*sqrt(6.) * x * (4.*z*z - x*x -y*y)
+        case(2) 
+          ylm =  ylm * 0.5*sqrt(15.) * z * (x*x -y*y)
+        case(3) 
+          ylm =  ylm * 0.25*sqrt(10.) * x * (x*x - 3. * y*y)
       end select
 
-!    case(2)
-!      select case (m)
-!        case(-2) 
-!          ylm =  ylm * sqrt(15.) * x * y * r_inv * r_inv
-!        case(-1) 
-!          ylm =  ylm * sqrt(15.) * y * z * r_inv * r_inv
-!        case(0) 
-!          ylm =  0.5 * ylm * sqrt(15.) * (2.*z*z - x*x - y*y) * r_inv * r_inv
-!        case(1) 
-!          ylm =  ylm * sqrt(15.) * z * x * r_inv * r_inv
-!        case(2) 
-!          ylm =  0.5 * ylm * sqrt(15.) * (x*x - y*y) * r_inv * r_inv
-!      end select
 
     case default
-      stop 'problem in Ylm of pseudo'
+      print *,  'l=', l
+      stop 'problem in Ylm of pseudo : Ylm not implemented (pseudo.irp.f)'
 
   end select