Break lines at 132 characters

Signed-off-by: Seyed Ali Ghasemi <info@gha3mi.com>
gha3mi · Jun 22, 2024 · 605394b · 605394b
1 parent dd7f080
commit 605394b
Show file tree

Hide file tree

Showing 8 changed files with 144 additions and 51 deletions.
diff --git a/example/nearest_point_1d.f90 b/example/nearest_point_1d.f90
@@ -48,7 +48,8 @@ program nearest_point_1d
     ! nearest_Xt: Corresponding parametric coordinates of the nearest point (optional)
     ! id: id of the nearest point (optional)
     call shape%nearest_point([4.5_rk, 4.5_rk, 5.0_rk], nearest_Xg, nearest_Xt, id)
-    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,a,1x,g0)','Nearest point on the curve:', nearest_Xg, ' with parametric coordinates:', nearest_Xt, ' and id:', id
+    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,a,1x,g0)',&
+        'Nearest point on the curve:', nearest_Xg, ' with parametric coordinates:', nearest_Xt, ' and id:', id
 
     !-----------------------------------------------------------------------------
     ! Nearest point on the curve (Optimization)
@@ -60,7 +61,8 @@ program nearest_point_1d
     ! nearest_Xt: Corresponding parametric coordinates of the nearest point
     ! nearest_Xg: Coordinates of the nearest point on the curve (optional)
     call shape%nearest_point2([4.5_rk, 4.5_rk, 5.0_rk], 1.0e-11_rk, 30, nearest_Xt, nearest_Xg)
-    print '(a,1x,g0,2x,g0,a,2x,g0,2x,g0)', 'Nearest point on the curve:', nearest_Xg, ' with parametric coordinates:', nearest_Xt
+    print '(a,1x,g0,2x,g0,a,2x,g0,2x,g0)',&
+        'Nearest point on the curve:', nearest_Xg, ' with parametric coordinates:', nearest_Xt
 
     !-----------------------------------------------------------------------------
     ! Finalizing

diff --git a/example/nearest_point_2d.f90 b/example/nearest_point_2d.f90
@@ -42,7 +42,8 @@ program nearest_point_2d
     ! nearest_Xt: Corresponding parametric coordinates of the nearest point (optional)
     ! id: id of the nearest point (optional)
     call shape%nearest_point([1.3_rk, 1.0_rk, 1.999999999_rk], nearest_Xg, nearest_Xt, id)
-    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0,2x,a,1x,g0)','Nearest point on the surface:', nearest_Xg, ' with parametric coordinates:', nearest_Xt, ' and id:', id
+    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0,2x,a,1x,g0)',&
+        'Nearest point on the surface:', nearest_Xg, ' with parametric coordinates:', nearest_Xt, ' and id:', id
 
     !-----------------------------------------------------------------------------
     ! Nearest point on the surface (Optimization)
@@ -54,7 +55,8 @@ program nearest_point_2d
     ! nearest_Xt: Corresponding parametric coordinates of the nearest point
     ! nearest_Xg: Coordinates of the nearest point on the surface (optional)
     call shape%nearest_point2([1.3_rk, 1.0_rk, 1.999999999_rk], 1.0e-11_rk, 30, nearest_Xt, nearest_Xg)
-    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0)', 'Nearest point on the surface:', nearest_Xg, ' with parametric coordinates:', nearest_Xt
+    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0)',&
+        'Nearest point on the surface:', nearest_Xg, ' with parametric coordinates:', nearest_Xt
 
     !-----------------------------------------------------------------------------
     ! Finalizing

diff --git a/example/nearest_point_3d.f90 b/example/nearest_point_3d.f90
@@ -27,7 +27,11 @@ program nearest_point_3d
     !> The weights of the control points (Wc) are optional.
     Wc = [1.0_rk, 1.1_rk, 1.11_rk, 1.0_rk, 0.5_rk, 0.5_rk, 1.2_rk, 1.0_rk]
 
-    call shape%set(knot1=[0.0_rk, 0.0_rk, 1.0_rk, 1.0_rk], knot2=[0.0_rk, 0.0_rk, 1.0_rk, 1.0_rk], knot3=[0.0_rk, 0.0_rk, 1.0_rk, 1.0_rk], Xc=Xc, Wc=Wc)
+    call shape%set(&
+        knot1=[0.0_rk, 0.0_rk, 1.0_rk, 1.0_rk],&
+        knot2=[0.0_rk, 0.0_rk, 1.0_rk, 1.0_rk],&
+        knot3=[0.0_rk, 0.0_rk, 1.0_rk, 1.0_rk],&
+        Xc=Xc, Wc=Wc)
 
     !-----------------------------------------------------------------------------
     ! Creating the NURBS volume
@@ -45,7 +49,8 @@ program nearest_point_3d
     ! nearest_Xt: Corresponding parametric coordinates of the nearest point (optional)
     ! id: id of the nearest point (optional)
     call shape%nearest_point([1.5_rk, 3.5_rk, 1.1_rk], nearest_Xg, nearest_Xt, id)
-    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0,2x,g0,2x,a,1x,g0)','Nearest point on the volume:', nearest_Xg, ' with parametric coordinates:', nearest_Xt, ' and id:', id
+    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0,2x,g0,2x,a,1x,g0)',&
+        'Nearest point on the volume:', nearest_Xg, ' with parametric coordinates:', nearest_Xt, ' and id:', id
 
     !-----------------------------------------------------------------------------
     ! Nearest point on the volume (Optimization)
@@ -57,7 +62,8 @@ program nearest_point_3d
     ! nearest_Xt: Corresponding parametric coordinates of the nearest point
     ! nearest_Xg: Coordinates of the nearest point on the volume (optional)
     call shape%nearest_point2([1.5_rk, 3.5_rk, 1.1_rk], 1.0e-11_rk, 500, nearest_Xt, nearest_Xg)
-    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0,2x,g0)', 'Nearest point on the volume:', nearest_Xg, ' with parametric coordinates:', nearest_Xt
+    print '(a,1x,g0,2x,g0,2x,g0,a,2x,g0,2x,g0,2x,g0)',&
+        'Nearest point on the volume:', nearest_Xg, ' with parametric coordinates:', nearest_Xt
 
     !-----------------------------------------------------------------------------
     ! Finalizing

diff --git a/src/forcad_nurbs_curve.f90 b/src/forcad_nurbs_curve.f90
@@ -1806,7 +1806,8 @@ impure subroutine nearest_point2(this, point_Xg, tol, maxit, nearest_Xt, nearest
 
             obj  = norm2(Xg - point_Xg) + 0.001_rk ! add a small number to avoid division by zero
             grad = dot_product((Xg-point_Xg)/obj, matmul(dTgc,this%Xc))
-            hess = dot_product(matmul(dTgc,this%Xc) - (Xg-point_Xg)/obj*grad, matmul(dTgc,this%Xc))/obj + dot_product((Xg-point_Xg)/obj, matmul(d2Tgc,this%Xc))
+            hess = dot_product(matmul(dTgc,this%Xc) - (Xg-point_Xg)/obj*grad, matmul(dTgc,this%Xc))/obj&
+                + dot_product((Xg-point_Xg)/obj, matmul(d2Tgc,this%Xc))
 
             ! debug
             print '(i3,1x,e20.10,1x,e20.10)', k, xk, abs(grad)

diff --git a/src/forcad_nurbs_surface.f90 b/src/forcad_nurbs_surface.f90
@@ -2431,10 +2431,18 @@ impure subroutine nearest_point2(this, point_Xg, tol, maxit, nearest_Xt, nearest
             grad(1) = dot_product((Xg-point_Xg)/obj, matmul(dTgc(:,1),this%Xc))
             grad(2) = dot_product((Xg-point_Xg)/obj, matmul(dTgc(:,2),this%Xc))
 
-            hess(1,1) = ( dot_product(matmul(dTgc(:,1),this%Xc), matmul(dTgc(:,1),this%Xc)) + dot_product((Xg-point_Xg), matmul(d2Tgc(1:this%nc(1)*this%nc(2)                        ,1),this%Xc)) ) /obj - ( dot_product(Xg-point_Xg, matmul(dTgc(:,1), this%Xc))*grad(1) ) / obj**2
-            hess(2,1) = ( dot_product(matmul(dTgc(:,1),this%Xc), matmul(dTgc(:,2),this%Xc)) + dot_product((Xg-point_Xg), matmul(d2Tgc(this%nc(1)*this%nc(2)+1:2*this%nc(1)*this%nc(2),1),this%Xc)) ) /obj - ( dot_product(Xg-point_Xg, matmul(dTgc(:,2), this%Xc))*grad(1) ) / obj**2
-            hess(1,2) = ( dot_product(matmul(dTgc(:,2),this%Xc), matmul(dTgc(:,1),this%Xc)) + dot_product((Xg-point_Xg), matmul(d2Tgc(1:this%nc(1)*this%nc(2)                        ,2),this%Xc)) ) /obj - ( dot_product(Xg-point_Xg, matmul(dTgc(:,1), this%Xc))*grad(2) ) / obj**2
-            hess(2,2) = ( dot_product(matmul(dTgc(:,2),this%Xc), matmul(dTgc(:,2),this%Xc)) + dot_product((Xg-point_Xg), matmul(d2Tgc(this%nc(1)*this%nc(2)+1:2*this%nc(1)*this%nc(2),2),this%Xc)) ) /obj - ( dot_product(Xg-point_Xg, matmul(dTgc(:,2), this%Xc))*grad(2) ) / obj**2
+            hess(1,1) = ( dot_product(matmul(dTgc(:,1),this%Xc), matmul(dTgc(:,1),this%Xc)) + &
+                dot_product((Xg-point_Xg), matmul(d2Tgc(1:this%nc(1)*this%nc(2)                        ,1),this%Xc)) ) /obj &
+                - ( dot_product(Xg-point_Xg, matmul(dTgc(:,1), this%Xc))*grad(1) ) / obj**2
+            hess(2,1) = ( dot_product(matmul(dTgc(:,1),this%Xc), matmul(dTgc(:,2),this%Xc)) + &
+                dot_product((Xg-point_Xg), matmul(d2Tgc(this%nc(1)*this%nc(2)+1:2*this%nc(1)*this%nc(2),1),this%Xc)) ) /obj &
+                - ( dot_product(Xg-point_Xg, matmul(dTgc(:,2), this%Xc))*grad(1) ) / obj**2
+            hess(1,2) = ( dot_product(matmul(dTgc(:,2),this%Xc), matmul(dTgc(:,1),this%Xc)) + &
+                dot_product((Xg-point_Xg), matmul(d2Tgc(1:this%nc(1)*this%nc(2)                        ,2),this%Xc)) ) /obj &
+                - ( dot_product(Xg-point_Xg, matmul(dTgc(:,1), this%Xc))*grad(2) ) / obj**2
+            hess(2,2) = ( dot_product(matmul(dTgc(:,2),this%Xc), matmul(dTgc(:,2),this%Xc)) + &
+                dot_product((Xg-point_Xg), matmul(d2Tgc(this%nc(1)*this%nc(2)+1:2*this%nc(1)*this%nc(2),2),this%Xc)) ) /obj &
+                - ( dot_product(Xg-point_Xg, matmul(dTgc(:,2), this%Xc))*grad(2) ) / obj**2
 
             ! debug
             print '(i3,1x,2e20.10,1x,e20.10)', k, xk, norm2(grad)
@@ -2776,10 +2784,18 @@ impure subroutine compute_d2Tgc_nurbs_2d_vector(Xt, knot1, knot2, degree, nc, ng
         d2Bi(1:nc(1)*nc(2)              ,2) = kron(dB2, dB1)
         d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2) = kron(d2B2, B1)
 
-        d2Tgc(i,1:nc(1)*nc(2)              ,1) = (d2Bi(1:nc(1)*nc(2)              ,1)*Wc - 2.0_rk*dTgc(i,:,1)*dot_product(dBi(:,1),Wc)                                 - Tgc(i,:)*dot_product(d2Bi(1:nc(1)*nc(2)              ,1),Wc)) / dot_product(Bi,Wc)
-        d2Tgc(i,nc(1)*nc(2)+1:2*nc(1)*nc(2),1) = (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1)*Wc - dTgc(i,:,1)*dot_product(dBi(:,2),Wc) - dTgc(i,:,2)*dot_product(dBi(:,1),Wc) - Tgc(i,:)*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1),Wc)) / dot_product(Bi,Wc)
-        d2Tgc(i,1:nc(1)*nc(2)              ,2) = (d2Bi(1:nc(1)*nc(2)              ,2)*Wc - dTgc(i,:,1)*dot_product(dBi(:,2),Wc) - dTgc(i,:,2)*dot_product(dBi(:,1),Wc) - Tgc(i,:)*dot_product(d2Bi(1:nc(1)*nc(2)              ,2),Wc)) / dot_product(Bi,Wc)
-        d2Tgc(i,nc(1)*nc(2)+1:2*nc(1)*nc(2),2) = (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2)*Wc - 2.0_rk*dTgc(i,:,2)*dot_product(dBi(:,2),Wc)                                 - Tgc(i,:)*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2),Wc)) / dot_product(Bi,Wc)
+        d2Tgc(i,1:nc(1)*nc(2)              ,1) = &
+            (d2Bi(1:nc(1)*nc(2)              ,1)*Wc - 2.0_rk*dTgc(i,:,1)*dot_product(dBi(:,1),Wc)                                 &
+            - Tgc(i,:)*dot_product(d2Bi(1:nc(1)*nc(2)              ,1),Wc)) / dot_product(Bi,Wc)
+        d2Tgc(i,nc(1)*nc(2)+1:2*nc(1)*nc(2),1) = &
+            (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1)*Wc - dTgc(i,:,1)*dot_product(dBi(:,2),Wc) - dTgc(i,:,2)*dot_product(dBi(:,1),Wc) &
+            - Tgc(i,:)*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1),Wc)) / dot_product(Bi,Wc)
+        d2Tgc(i,1:nc(1)*nc(2)              ,2) = &
+            (d2Bi(1:nc(1)*nc(2)              ,2)*Wc - dTgc(i,:,1)*dot_product(dBi(:,2),Wc) - dTgc(i,:,2)*dot_product(dBi(:,1),Wc) &
+            - Tgc(i,:)*dot_product(d2Bi(1:nc(1)*nc(2)              ,2),Wc)) / dot_product(Bi,Wc)
+        d2Tgc(i,nc(1)*nc(2)+1:2*nc(1)*nc(2),2) = &
+            (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2)*Wc - 2.0_rk*dTgc(i,:,2)*dot_product(dBi(:,2),Wc)                                 &
+            - Tgc(i,:)*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2),Wc)) / dot_product(Bi,Wc)
     end do
 end subroutine
 !===============================================================================
@@ -2829,10 +2845,18 @@ impure subroutine compute_d2Tgc_nurbs_2d_scalar(Xt, knot1, knot2, degree, nc, Wc
     d2Bi(1:nc(1)*nc(2)              ,2) = kron(dB2, dB1)
     d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2) = kron(d2B2, B1)
 
-    d2Tgc(1:nc(1)*nc(2)              ,1) = (d2Bi(1:nc(1)*nc(2)              ,1)*Wc - 2.0_rk*dTgc(:,1)*dot_product(dBi(:,1),Wc)                               - Tgc*dot_product(d2Bi(1:nc(1)*nc(2)              ,1),Wc)) / dot_product(Bi,Wc)
-    d2Tgc(nc(1)*nc(2)+1:2*nc(1)*nc(2),1) = (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1)*Wc - dTgc(:,1)*dot_product(dBi(:,2),Wc) - dTgc(:,2)*dot_product(dBi(:,1),Wc) - Tgc*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1),Wc)) / dot_product(Bi,Wc)
-    d2Tgc(1:nc(1)*nc(2)              ,2) = (d2Bi(1:nc(1)*nc(2)              ,2)*Wc - dTgc(:,1)*dot_product(dBi(:,2),Wc) - dTgc(:,2)*dot_product(dBi(:,1),Wc) - Tgc*dot_product(d2Bi(1:nc(1)*nc(2)              ,2),Wc)) / dot_product(Bi,Wc)
-    d2Tgc(nc(1)*nc(2)+1:2*nc(1)*nc(2),2) = (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2)*Wc - 2.0_rk*dTgc(:,2)*dot_product(dBi(:,2),Wc)                               - Tgc*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2),Wc)) / dot_product(Bi,Wc)
+    d2Tgc(1:nc(1)*nc(2)              ,1) = &
+        (d2Bi(1:nc(1)*nc(2)              ,1)*Wc - 2.0_rk*dTgc(:,1)*dot_product(dBi(:,1),Wc)                               &
+        - Tgc*dot_product(d2Bi(1:nc(1)*nc(2)              ,1),Wc)) / dot_product(Bi,Wc)
+    d2Tgc(nc(1)*nc(2)+1:2*nc(1)*nc(2),1) = &
+        (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1)*Wc - dTgc(:,1)*dot_product(dBi(:,2),Wc) - dTgc(:,2)*dot_product(dBi(:,1),Wc) &
+        - Tgc*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),1),Wc)) / dot_product(Bi,Wc)
+    d2Tgc(1:nc(1)*nc(2)              ,2) = &
+        (d2Bi(1:nc(1)*nc(2)              ,2)*Wc - dTgc(:,1)*dot_product(dBi(:,2),Wc) - dTgc(:,2)*dot_product(dBi(:,1),Wc) &
+        - Tgc*dot_product(d2Bi(1:nc(1)*nc(2)              ,2),Wc)) / dot_product(Bi,Wc)
+    d2Tgc(nc(1)*nc(2)+1:2*nc(1)*nc(2),2) = &
+        (d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2)*Wc - 2.0_rk*dTgc(:,2)*dot_product(dBi(:,2),Wc)                               &
+        - Tgc*dot_product(d2Bi(nc(1)*nc(2)+1:2*nc(1)*nc(2),2),Wc)) / dot_product(Bi,Wc)
 end subroutine
 !===============================================================================