赫斯史密斯面元法，求解无界流中物体的附加质量

crazybiscuit

已于 2025-02-24 10:20:24 修改

阅读量608

点赞数 6

文章标签：开源软件

于 2024-06-16 21:40:07 首次发布

本文链接：https://blog.youkuaiyun.com/crazybiscuit/article/details/139726416

版权

一、代码

本代码使用Fortran语言编写。在理论上，使用四边形网格、使用源分布模型，用于计算无界流中单浮体或多浮体情况下的目标浮体的附加质量。

详细理论参考《船舶在波浪中运动的势流理论》（戴遗山段文洋著）及相关MOOC。

program Hess_Smith_MultiFloats
implicit none

interface
  subroutine wenjian(name,Npoints,Npanels,NpointsT,NpanelsT,coor,ele,AVO,nj,n_xi,n_eta,coor_local,area)
    character*12 :: name
    integer :: Npoints,Npanels
    integer :: NpointsT,NpanelsT
    real*8,allocatable  :: coor(:,:)
    integer,allocatable :: ele(:,:)
    real*8,allocatable  :: AVO(:,:)                !面元原点平均坐标
    real*8,allocatable  :: nj(:,:)                 !面元单位法向矢量(大地系)
    real*8,allocatable  :: n_xi(:,:)               !面元局部坐标系xi方向的单位向量（大地系）
    real*8,allocatable  :: n_eta(:,:)              !面元局部坐标系eta方向的单位向量（大地系）
    real*8,allocatable  :: coor_local(:,:)         !投影后点的局部坐标
    real*8,allocatable  :: area(:)                 !每个面元的面积
    real*8  :: brick=0.0                           !哪里需要就用在哪
    real*8  :: brick1=0.0                          !哪里需要就用在哪
  end subroutine
end interface

interface
  subroutine yingxiangxishu(Npoints,Npanels,AVO,nj,n_xi,n_eta,aij,cij,coor_local)
    integer :: Npoints,Npanels
    real*8,allocatable :: AVO(:,:)                 !面元原点平均坐标
    real*8,allocatable :: nj(:,:)                  !面元单位法向矢量(大地系)
    real*8,allocatable :: n_xi(:,:)                !面元局部坐标系xi方向的单位向量（大地系）
    real*8,allocatable :: n_eta(:,:)               !面元局部坐标系eta方向的单位向量（大地系）
    real*8,allocatable :: aij(:,:)                 !影响系数矩阵
    real*8,allocatable :: cij(:,:)
    real*8,allocatable :: coor_local(:,:)          !投影后点的局部坐标
  end subroutine

  subroutine AdditionalQuality(Npoints,Npanels,NpointsT,NpanelsT,AVO,nj,aij,cij,area,ni,bi,sigma,phi_pi,mji)
    integer :: Npoints,Npanels
    integer :: NpointsT,NpanelsT
    real*8,allocatable :: AVO(:,:)          !面元原点平均坐标
    real*8,allocatable :: nj(:,:)           !面元单位法向矢量(大地系)
    real*8,allocatable :: aij(:,:)          !影响系数矩阵
    real*8,allocatable :: cij(:,:)
    real*8,allocatable :: area(:)           !每个面元的面积
    real*8,allocatable :: ni(:,:)
    real*8,allocatable :: bi(:)
    real*8,allocatable :: sigma(:)
    real*8,allocatable :: phi_pi(:)
    real*8 :: mji(6,6)
  end subroutine
end interface


character*12 :: name
integer :: i,j
integer :: Npoints,Npanels             !总网格数，如果是单浮体那就=NpointT
integer :: NpointsT,NpanelsT           !NpointsTarget,目标浮体的网格数
real*8,allocatable :: coor(:,:)
integer,allocatable :: ele(:,:)
real*8,allocatable :: AVO(:,:)         !面元原点平均坐标
real*8,allocatable :: nj(:,:)          !面元单位法向矢量(大地系)
real*8,allocatable :: n_xi(:,:)        !面元局部坐标系xi方向的单位向量（大地系）
real*8,allocatable :: n_eta(:,:)       !面元局部坐标系eta方向的单位向量（大地系）
real*8,allocatable :: aij(:,:)         !影响系数矩阵
real*8,allocatable :: cij(:,:)
real*8,allocatable :: coor_local(:,:)  !投影后点的局部坐标
real*8,allocatable :: area(:)          !每个面元的面积
real*8,allocatable :: ni(:,:)
real*8,allocatable :: bi(:)
real*8,allocatable :: sigma(:)
real*8,allocatable :: phi_pi(:)
real*8 :: mji(6,6)


print *, '程序正在运行，请稍等……'


call wenjian(name,Npoints,Npanels,NpointsT,NpanelsT,coor,ele,AVO,nj,n_xi,n_eta,coor_local,area)
call yingxiangxishu(Npoints,Npanels,AVO,nj,n_xi,n_eta,aij,cij,coor_local)
call AdditionalQuality(Npoints,Npanels,NpointsT,NpanelsT,AVO,nj,aij,cij,area,ni,bi,sigma,phi_pi,mji)


open(20,file='mji.txt')
do j=1,6
  do i=1,5
   write(20,"(F12.3)",advance='no') mji(j,i)
  end do
  write(20,"(F12.3)") mji(j,6)
end do
close(20)

print *, 'm11=',mji(1,1)
print *, 'm33=',mji(3,3)


print *, '计算结束，请关闭'
read(*,*) i

end program


!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!读取文件并处理相关数据
subroutine wenjian(name,Npoints,Npanels,NpointsT,NpanelsT,coor,ele,AVO,nj,n_xi,n_eta,coor_local,area)
implicit none

character*12 :: name
integer :: i,j
integer :: Npoints,Npanels
integer :: NpointsT,NpanelsT
real*8,allocatable  :: coor(:,:)
integer,allocatable :: ele(:,:)
real*8,allocatable  :: AVO(:,:)                !面元原点平均坐标
real*8,allocatable  :: nj(:,:)                 !面元单位法向矢量(大地系)
real*8,allocatable  :: n_xi(:,:)               !面元局部坐标系xi方向的单位向量（大地系）
real*8,allocatable  :: n_eta(:,:)              !面元局部坐标系eta方向的单位向量（大地系）
real*8  :: V13(3)=0.0
real*8  :: V24(3)=0.0
real*8  :: Oq(3)=(/0.0,0.0,0.0/)
real*8  :: Oq1(3)=(/0.0,0.0,0.0/)
real*8,allocatable  :: coor_local(:,:)            !投影后点的局部坐标
real*8,allocatable  :: area(:)                    !每个面元的面积
real*8  :: brick=0.0                  !哪里需要就用在哪
real*8  :: brick1=0.0                 !哪里需要就用在哪


open(10,file="input.txt")
read(10,*) Npoints,Npanels,NpointsT,NpanelsT
allocate(coor(Npoints,3))
allocate(ele(Npanels,4))
do i=1,Npoints
  read(10,*) coor(i,1),coor(i,2),coor(i,3)
end do
do i=1,Npanels
  read(10,*) ele(i,1),ele(i,2),ele(i,3),ele(i,4)           !为了满足求Sx,Sy,Sz时逆时针排列（平面格林公式，法向和线积分方向的关系）
end do
close(10)


!求面元原点平均坐标
allocate(AVO(Npanels,3))

do i=1,Npanels
  AVO(i,1)=(coor(ele(i,1),1)+coor(ele(i,2),1)+coor(ele(i,3),1)+coor(ele(i,4),1))/4.0
  AVO(i,2)=(coor(ele(i,1),2)+coor(ele(i,2),2)+coor(ele(i,3),2)+coor(ele(i,4),2))/4.0
  AVO(i,3)=(coor(ele(i,1),3)+coor(ele(i,2),3)+coor(ele(i,3),3)+coor(ele(i,4),3))/4.0
end do

!求面元单位法向矢量(大地系)
allocate(nj(Npanels,3))
do i=1,Npanels
  V13(1)=coor(ele(i,3),1)-coor(ele(i,1),1)  !每个ele的1、3点形成的向量
  V13(2)=coor(ele(i,3),2)-coor(ele(i,1),2)
  V13(3)=coor(ele(i,3),3)-coor(ele(i,1),3)
  V24(1)=coor(ele(i,4),1)-coor(ele(i,2),1)  !每个ele的2、4点形成的向量
  V24(2)=coor(ele(i,4),2)-coor(ele(i,2),2)
  V24(3)=coor(ele(i,4),3)-coor(ele(i,2),3)
  brick=SQRT((V13(2)*V24(3)-V24(2)*V13(3))**2.0+(V24(1)*V13(3)-V13(1)*V24(3))**2.0+(V13(1)*V24(2)-V24(1)*V13(2))**2.0)
  
  nj(i,1)=(V13(2)*V24(3)-V24(2)*V13(3))/brick  !brick是向量13与24叉乘结果的模
  nj(i,2)=(V24(1)*V13(3)-V13(1)*V24(3))/brick
  nj(i,3)=(V13(1)*V24(2)-V24(1)*V13(2))/brick

end do


!求每个面元4个点投影的局部系坐标（xi（即kesei）和eta）
allocate(n_xi(Npanels,3))
allocate(n_eta(Npanels,3))
allocate(coor_local(4*Npanels,3))
allocate(area(Npanels))
do i=1,Npanels
  do j=1,4
    Oq(1)=coor(ele(i,j),1)-AVO(i,1)
    Oq(2)=coor(ele(i,j),2)-AVO(i,2)
    Oq(3)=coor(ele(i,j),3)-AVO(i,3)
    brick=Oq(1)*nj(i,1)+Oq(2)*nj(i,2)+Oq(3)*nj(i,3)  !brick:Oq与法向的点乘

    Oq1(1)=Oq(1)-brick*nj(i,1)
    Oq1(2)=Oq(2)-brick*nj(i,2)
    Oq1(3)=Oq(3)-brick*nj(i,3)


    if(j==1) then      !求面元局部坐标系xi方向的单位向量（大地系）(认为每个ele的第1个点是xi轴方向)
      brick1=SQRT(Oq1(1)**2.0+Oq1(2)**2.0+Oq1(3)**2.0)  !brick1是向量Oq1的模
      n_xi(i,1)=Oq1(1)/brick1
      n_xi(i,2)=Oq1(2)/brick1
      n_xi(i,3)=Oq1(3)/brick1

      n_eta(i,1)=nj(i,2)*n_xi(i,3)-nj(i,3)*n_xi(i,2)   !叉乘得到eta（大地系）
      n_eta(i,2)=nj(i,3)*n_xi(i,1)-nj(i,1)*n_xi(i,3)
      n_eta(i,3)=nj(i,1)*n_xi(i,2)-nj(i,2)*n_xi(i,1)
    end if


    coor_local(4*(i-1)+j,1)=Oq1(1)*n_xi(i,1)+Oq1(2)*n_xi(i,2)+Oq1(3)*n_xi(i,3)   !求每个ele的4个点投影的局部坐标
    coor_local(4*(i-1)+j,2)=Oq1(1)*n_eta(i,1)+Oq1(2)*n_eta(i,2)+Oq1(3)*n_eta(i,3)
    coor_local(4*(i-1)+j,3)=Oq1(1)*nj(i,1)+Oq1(2)*nj(i,2)+Oq1(3)*nj(i,3)         !!!应该不用求，就是0才对

  end do


  !求每个面元的面积（分成两个三角形）
  area(i)=0.5*ABS((coor_local(4*(i-1)+4,1)-coor_local(4*(i-1)+1,1))*(coor_local(4*(i-1)+2,2)-coor_local(4*(i-1)+1,2))-&
  (coor_local(4*(i-1)+4,2)-coor_local(4*(i-1)+1,2))*(coor_local(4*(i-1)+2,1)-coor_local(4*(i-1)+1,1)))+&
  0.5*ABS((coor_local(4*(i-1)+4,1)-coor_local(4*(i-1)+3,1))*(coor_local(4*(i-1)+2,2)-coor_local(4*(i-1)+3,2))&
  -(coor_local(4*(i-1)+4,2)-coor_local(4*(i-1)+3,2))*(coor_local(4*(i-1)+2,1)-coor_local(4*(i-1)+3,1)))

end do

end subroutine



!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!计算Sx,Sy,Sz,S,aij,cij
subroutine yingxiangxishu(Npoints,Npanels,AVO,nj,n_xi,n_eta,aij,cij,coor_local)
implicit none

real*8,parameter :: pi=4.0*ATAN(1.0)      !π
integer :: i,j,k
integer :: Npoints,Npanels
real*8,allocatable :: AVO(:,:)  !面元原点平均坐标
real*8,allocatable :: nj(:,:)   !面元单位法向矢量(大地系)
real*8,allocatable :: n_xi(:,:) !面元局部坐标系xi方向的单位向量（大地系）
real*8,allocatable :: n_eta(:,:)!面元局部坐标系eta方向的单位向量（大地系）
real*8,allocatable :: aij(:,:)  !影响系数矩阵
real*8,allocatable :: cij(:,:)
real*8 :: AVO_local(3)=(/0.0,0.0,0.0/)
real*8 :: Oq1(3)=(/0.0,0.0,0.0/)
real*8,allocatable :: coor_local(:,:)  !投影后点的局部坐标
real*8 :: l(4),r(4),c(4),h(4),m(4)     !求Sx,Sy,Sz,S用
real*8 :: Sx=0.0,Sy=0.0,Sz=0.0,S=0.0
real*8 :: Sx1=0.0,Sy1=0.0,Sz1=0.0



allocate(aij(Npanels,Npanels))
allocate(cij(Npanels,Npanels))
do j=1,Npanels
  do i=1,Npanels
      Oq1(1)=AVO(i,1)-AVO(j,1)                !求第i个面元坐标原点在第j个面元的局部系坐标AVO_local
      Oq1(2)=AVO(i,2)-AVO(j,2)
      Oq1(3)=AVO(i,3)-AVO(j,3)

      AVO_local(1)=Oq1(1)*n_xi(j,1)+Oq1(2)*n_xi(j,2)+Oq1(3)*n_xi(j,3)
      AVO_local(2)=Oq1(1)*n_eta(j,1)+Oq1(2)*n_eta(j,2)+Oq1(3)*n_eta(j,3)
      AVO_local(3)=Oq1(1)*nj(j,1)+Oq1(2)*nj(j,2)+Oq1(3)*nj(j,3)           !到这里求完AVO_local

      do k=1,3         !求l
        l(k)=SQRT((coor_local(4*(j-1)+k+1,1)-coor_local(4*(j-1)+k,1))**2.0+(coor_local(4*(j-1)+k+1,2)-coor_local(4*(j-1)+k,2))**2.0)
      end do
      l(4)=SQRT((coor_local(4*(j-1)+1,1)-coor_local(4*(j-1)+4,1))**2.0+(coor_local(4*(j-1)+1,2)-coor_local(4*(j-1)+4,2))**2.0)


      do k=1,4
        c(k)=(coor_local(4*(j-1)+k,1)-AVO_local(1))**2.0+AVO_local(3)**2.0    !求c
        r(k)=SQRT(c(k)+(coor_local(4*(j-1)+k,2)-AVO_local(2))**2.0)           !求r
        h(k)=(coor_local(4*(j-1)+k,1)-AVO_local(1))*(coor_local(4*(j-1)+k,2)-AVO_local(2))  !求h
      end do
      
      do k=1,3    !求m
        m(k)=(coor_local(4*(j-1)+k+1,2)-coor_local(4*(j-1)+k,2))/(coor_local(4*(j-1)+k+1,1)-coor_local(4*(j-1)+k,1))
      end do
      m(4)=(coor_local(4*(j-1)+1,2)-coor_local(4*(j-1)+4,2))/(coor_local(4*(j-1)+1,1)-coor_local(4*(j-1)+4,1))




      Sx=0.0      !求Sx
      do k=1,3
        if(ABS((coor_local(4*(j-1)+k+1,2)-coor_local(4*(j-1)+k,2)))<=1.0E-12) then
            Sx=Sx
        else
          Sx=Sx-(coor_local(4*(j-1)+k+1,2)-coor_local(4*(j-1)+k,2))/l(k)*LOG((r(k)+r(k+1)+l(k))/(r(k)+r(k+1)-l(k)))
        end if
      end do
      if(ABS((coor_local(4*(j-1)+4,2)-coor_local(4*(j-1)+1,2)))<=1.0E-12) then
        Sx=Sx
      else
        Sx=Sx-(coor_local(4*(j-1)+1,2)-coor_local(4*(j-1)+4,2))/l(4)*LOG((r(4)+r(1)+l(4))/(r(4)+r(1)-l(4)))
      end if


      Sy=0.0      !求Sy
      do k=1,3
        if(ABS((coor_local(4*(j-1)+k+1,1)-coor_local(4*(j-1)+k,1)))<=1.0E-12) then
          Sy=Sy
        else
          Sy=Sy+(coor_local(4*(j-1)+k+1,1)-coor_local(4*(j-1)+k,1))/l(k)*LOG((r(k)+r(k+1)+l(k))/(r(k)+r(k+1)-l(k)))
        end if
      end do
      if(ABS((coor_local(4*(j-1)+4,1)-coor_local(4*(j-1)+1,1)))<=1.0E-12) then
        Sy=Sy
      else
        Sy=Sy+(coor_local(4*(j-1)+1,1)-coor_local(4*(j-1)+4,1))/l(4)*LOG((r(4)+r(1)+l(4))/(r(4)+r(1)-l(4)))
      end if


      Sz=0.0      !求Sz
      if(ABS(AVO_local(3))<=1.0E-12) then
        Sz=0.0
      else
        do k=1,3
          if(ABS((coor_local(4*(j-1)+k+1,1)-coor_local(4*(j-1)+k,1)))<=1.0E-12) then    !书上P27页xi_i+1=xi_i的情况：对应线段上Sz=0
            Sz=Sz
          else
            Sz=Sz+(ATAN((m(k)*c(k)-h(k))/AVO_local(3)/r(k))-ATAN((m(k)*c(k+1)-h(k+1))/AVO_local(3)/r(k+1)))
          end if
        end do
        if(ABS((coor_local(4*(j-1)+4,1)-coor_local(4*(j-1)+1,1)))<=1.0E-12) then
          Sz=Sz
        else
          Sz=Sz+(ATAN((m(4)*c(4)-h(4))/AVO_local(3)/r(4))-ATAN((m(4)*c(1)-h(1))/AVO_local(3)/r(1)))
        end if
      end if

      S=0.0      !求S
      do k=1,3
        S=S+((coor_local(4*(j-1)+k+1,1)-coor_local(4*(j-1)+k,1))*(AVO_local(2)-coor_local(4*(j-1)+k,2))&
        -(coor_local(4*(j-1)+k+1,2)-coor_local(4*(j-1)+k,2))*(AVO_local(1)-coor_local(4*(j-1)+k,1)))&
        /l(k)*LOG((r(k)+r(k+1)+l(k))/(r(k)+r(k+1)-l(k)))
      end do
      S=S+((coor_local(4*(j-1)+1,1)-coor_local(4*(j-1)+4,1))*(AVO_local(2)-coor_local(4*(j-1)+4,2))&
      -(coor_local(4*(j-1)+1,2)-coor_local(4*(j-1)+4,2))*(AVO_local(1)-coor_local(4*(j-1)+4,1)))&
      /l(4)*LOG((r(4)+r(1)+l(4))/(r(4)+r(1)-l(4)))
      S=S+AVO_local(3)*Sz

      cij(i,j)=S    !求cij
    

    if(i==j) then
      aij(i,j)=2.0*pi
    else if(i/=j) then                                  !将（Sx,Sy,Sz）转回大地坐标系
      Sx1=Sx*n_xi(j,1)+Sy*n_eta(j,1)+Sz*nj(j,1)
      Sy1=Sx*n_xi(j,2)+Sy*n_eta(j,2)+Sz*nj(j,2)
      Sz1=Sx*n_xi(j,3)+Sy*n_eta(j,3)+Sz*nj(j,3)

      aij(i,j)=Sx1*nj(j,1)+Sy1*nj(j,2)+Sz1*nj(j,3)     !求aij
    end if
  end do
end do

end subroutine




!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!计算附加质量
subroutine AdditionalQuality(Npoints,Npanels,NpointsT,NpanelsT,AVO,nj,aij,cij,area,ni,bi,sigma,phi_pi,mji)
implicit none

real*8,parameter :: rho=1000.0        !水的密度
integer :: i,j,k,k1
integer :: Npoints,Npanels
integer :: NpointsT,NpanelsT
real*8,allocatable :: AVO(:,:)  !面元原点平均坐标
real*8,allocatable :: nj(:,:)   !面元单位法向矢量(大地系)
real*8,allocatable :: aij(:,:)  !影响系数矩阵
real*8,allocatable :: cij(:,:)
real*8,allocatable :: area(:)        !每个面元的面积
real*8,allocatable :: ni(:,:)
real*8,allocatable :: bi(:)
real*8,allocatable :: sigma(:)
real*8,allocatable :: phi_pi(:)
real*8 :: mji(6,6)



!计算bi，并求解线代方程组aij*sigma=bi(bi实际上是ni不同辐射势1-6)
allocate(ni(Npanels,6))
do i=1,Npanels
    ni(i,1)=nj(i,1)                        !根据式（2.1.17）、（2.1.15)
    ni(i,2)=nj(i,2)
    ni(i,3)=nj(i,3)
    ni(i,4)=AVO(i,2)*nj(i,3)-AVO(i,3)*nj(i,2)
    ni(i,5)=AVO(i,3)*nj(i,1)-AVO(i,1)*nj(i,3)
    ni(i,6)=AVO(i,1)*nj(i,2)-AVO(i,2)*nj(i,1)
end do

mji=0.0
allocate(sigma(Npanels))
allocate(bi(Npanels))
allocate(phi_pi(Npanels))
do k=1,6                        !求k方向自由度的运动在k1方向产生的附加质量mji(i对应k，j对应k1)
  do i=1,Npanels
    bi(i)=ni(i,k)
  end do

  sigma=0.0
  call gaussRK(aij,bi,sigma,Npanels)

  phi_pi(1:Npanels)=0.0         !式(2.7.20)求的：k方向自由度运动时，所有面元的分布源在pi点的诱导速度势
  do i=1,NpanelsT
    do j=1,Npanels
      phi_pi(i)=phi_pi(i)+cij(i,j)*sigma(j)
    end do
  end do

  do k1=1,6                                                  !求mji
    do i=1,NpanelsT
      mji(k1,k)=mji(k1,k)+rho*phi_pi(i)*ni(i,k1)*area(i)     !i对应k，j对应k1
    end do
  end do

end do


end subroutine


!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!高斯消去法（doolittle分解）
subroutine gaussRK(G,b,a,n)
implicit real*8(a-z)

integer :: i,j,k
integer :: n           !格拉姆是n阶方阵
real*8  :: G(n,n)      !格拉姆矩阵   Ax=b中的A
real*8  :: b(n)        !Ax=b中的b
real*8  :: a(n)        !拟合的函数的基phi的系数     Ax=b中的x

real*8,allocatable  :: m(:)    !见下方解释    m=Ua
real*8,allocatable  :: L(:,:)     
real*8,allocatable  :: U(:,:)
real*8 :: sigma=0.             !用于累加求和

allocate(m(n))
allocate(L(n,n))
allocate(U(n,n))

!本题中doolittle分解：Ga=b,G=LU,LUa=b,Ua=m,Lm=b


!求L、U、m、a(见P153页)
do i=1,n
  L(i,i)=1.0
end do

do i=1,n
  U(1,i)=G(1,i)
end do

do i=2,n
  L(i,1)=G(i,1)/U(1,1)
end do

do i=2,n
  do j=i,n
    do k=1,i-1
      sigma=sigma+L(i,k)*U(k,j)
    end do
    U(i,j)=G(i,j)-sigma
    sigma=0.
  end do
  
  if(i==n)  cycle
  do j=i+1,n
    do k=1,i-1
      sigma=sigma+L(j,k)*U(k,i)
    end do
    L(j,i)=(G(j,i)-sigma)/U(i,i)
    sigma=0.
  end do
end do

m(1)=b(1)
do i=2,n
  do k=1,i-1
    sigma=sigma+L(i,k)*m(k)
  end do
  m(i)=b(i)-sigma
  sigma=0.
end do

a(n)=m(n)/U(n,n)
do i=n-1,1,-1
  do k=i+1,n
    sigma=sigma+U(i,k)*a(k)
  end do
  a(i)=(m(i)-sigma)/U(i,i)
  sigma=0.
end do

end subroutine

二、输入文件

作为输入的“input.txt”文件有下述要求。

1、文件的内容参考代码中读取文件部分。

2、每个网格（ele）有4个点，input文件中点的排序应保证：按1234的顺序用右手螺旋定则判断方向，方向应指向物体外侧。（在代码中有使用点13和点24的向量叉乘给定面元法向的部分，但这部分代码好像不对，和理论正好写反了，我需要进一步研究。但是计算结果是正确的。）

3、如果计算单浮体，第一行仍要有4个数字，只是前两个重复、后两个也重复。

4、如果计算多浮体，请确保目标浮体的点（points）和网格（panels或者说ele）都在各自位置的最前面。以下面这个双浮体的input文件为例：372个Npoints节点中，前186个为目标浮体的节点NpointsT，剩余的为另一个浮体的节点；368个Npanels中，前184个为目标浮体的网格NpanelsT.

372 368   186 184      !分别是Npoints,Npanels,NpointsT,NpanelsT
  -4.7182400E+00   -1.8171000E-01   -2.7587000E-01
  -4.8324000E+00   -2.1791000E-01   -1.3313000E-01
  -4.9673800E+00    0.0000000E+00   -1.1169000E-01
  -4.8267600E+00    0.0000000E+00   -2.5997000E-01
  -4.2278100E+00    0.0000000E+00   -5.3229000E-01
  -4.2410100E+00   -2.8370000E-01   -4.4400000E-01
  -4.2426500E+00   -4.4967000E-01   -2.7299000E-01
  -3.5547400E+00   -5.9588000E-01   -3.6771000E-01
  -3.5400200E+00   -7.0611000E-01    0.0000000E+00
  -4.2279000E+00   -5.3226000E-01    0.0000000E+00
  -4.8387500E+00   -2.5076000E-01    0.0000000E+00
  -3.5550600E+00   -3.7059000E-01   -5.9401000E-01
  -2.8492200E+00   -4.3026000E-01   -6.9701000E-01
  -2.8491700E+00   -6.9697000E-01   -4.3021000E-01
  -3.5407600E+00    0.0000000E+00   -7.0595000E-01
  -2.8438200E+00    0.0000000E+00   -8.2126000E-01
  -2.1367200E+00    0.0000000E+00   -9.0282000E-01
  -2.1404800E+00   -4.7009000E-01   -7.6877000E-01
  -1.4284700E+00   -4.9948000E-01   -8.1479000E-01
  -1.4282600E+00   -8.1543000E-01   -4.9756000E-01
  -2.1403300E+00   -7.6784000E-01   -4.7067000E-01
  -1.4260000E+00    0.0000000E+00   -9.5798000E-01
  -7.1305000E-01    0.0000000E+00   -9.8871000E-01
  -7.1466000E-01   -5.1691000E-01   -8.4100000E-01
   0.0000000E+00   -5.1135000E-01   -8.5791000E-01
   0.0000000E+00   -8.5971000E-01   -5.0827000E-01
  -7.1469000E-01   -8.4232000E-01   -5.1482000E-01
   0.0000000E+00    0.0000000E+00   -1.0000000E+00
   0.0000000E+00   -1.0000000E+00    0.0000000E+00
  -7.1311000E-01   -9.8871000E-01    0.0000000E+00
  -1.4258100E+00   -9.5800000E-01    0.0000000E+00
  -2.1366800E+00   -9.0283000E-01    0.0000000E+00
  -2.8440800E+00   -8.2123000E-01    0.0000000E+00
  -5.0000000E+00    0.0000000E+00    0.0000000E+00
   4.8267600E+00    0.0000000E+00   -2.5997000E-01
   4.9673800E+00    0.0000000E+00   -1.1169000E-01
   4.8324000E+00   -2.1791000E-01   -1.3313000E-01
   4.7182400E+00   -1.8171000E-01   -2.7587000E-01
   4.2410100E+00   -2.8370000E-01   -4.4400000E-01
   4.2278100E+00    0.0000000E+00   -5.3229000E-01
   4.2426500E+00   -4.4967000E-01   -2.7299000E-01
   4.2279000E+00   -5.3226000E-01    0.0000000E+00
   3.5400200E+00   -7.0611000E-01    0.0000000E+00
   3.5547400E+00   -5.9588000E-01   -3.6771000E-01
   4.8387500E+00   -2.5076000E-01    0.0000000E+00
   3.5550600E+00   -3.7059000E-01   -5.9401000E-01
   2.8491700E+00   -6.9697000E-01   -4.3021000E-01
   2.8492200E+00   -4.3026000E-01   -6.9701000E-01
   3.5407600E+00    0.0000000E+00   -7.0595000E-01
   2.8438200E+00    0.0000000E+00   -8.2126000E-01
   2.1404800E+00   -4.7009000E-01   -7.6877000E-01
   2.1367200E+00    0.0000000E+00   -9.0282000E-01
   2.1403300E+00   -7.6784000E-01   -4.7067000E-01
   1.4282600E+00   -8.1543000E-01   -4.9756000E-01
   1.4284700E+00   -4.9948000E-01   -8.1479000E-01
   1.4260000E+00    0.0000000E+00   -9.5798000E-01
   7.1466000E-01   -5.1691000E-01   -8.4100000E-01
   7.1305000E-01    0.0000000E+00   -9.8871000E-01
   7.1469000E-01   -8.4232000E-01   -5.1482000E-01
   7.1311000E-01   -9.8871000E-01    0.0000000E+00
   1.4258100E+00   -9.5800000E-01    0.0000000E+00
   2.1366800E+00   -9.0283000E-01    0.0000000E+00
   2.8440800E+00   -8.2123000E-01    0.0000000E+00
   5.0000000E+00    0.0000000E+00    0.0000000E+00
  -4.8324000E+00    2.1791000E-01   -1.3313000E-01
  -4.7182400E+00    1.8171000E-01   -2.7587000E-01
  -4.2410100E+00    2.8370000E-01   -4.4400000E-01
  -4.2426500E+00    4.4967000E-01   -2.7299000E-01
  -4.2279000E+00    5.3226000E-01    0.0000000E+00
  -3.5400200E+00    7.0611000E-01    0.0000000E+00
  -3.5547400E+00    5.9588000E-01   -3.6771000E-01
  -4.8387500E+00    2.5076000E-01    0.0000000E+00
  -3.5550600E+00    3.7059000E-01   -5.9401000E-01
  -2.8491700E+00    6.9697000E-01   -4.3021000E-01
  -2.8492200E+00    4.3026000E-01   -6.9701000E-01
  -2.1404800E+00    4.7009000E-01   -7.6877000E-01
  -2.1403300E+00    7.6784000E-01   -4.7067000E-01
  -1.4282600E+00    8.1543000E-01   -4.9756000E-01
  -1.4284700E+00    4.9948000E-01   -8.1479000E-01
  -7.1466000E-01    5.1691000E-01   -8.4100000E-01
  -7.1469000E-01    8.4232000E-01   -5.1482000E-01
   0.0000000E+00    8.5971000E-01   -5.0827000E-01
   0.0000000E+00    5.1135000E-01   -8.5791000E-01
  -7.1311000E-01    9.8871000E-01    0.0000000E+00
   0.0000000E+00    1.0000000E+00    0.0000000E+00
  -1.4258100E+00    9.5800000E-01    0.0000000E+00
  -2.1366800E+00    9.0283000E-01    0.0000000E+00
  -2.8440800E+00    8.2123000E-01    0.0000000E+00
   4.7182400E+00    1.8171000E-01   -2.7587000E-01
   4.8324000E+00    2.1791000E-01   -1.3313000E-01
   4.2410100E+00    2.8370000E-01   -4.4400000E-01
   4.2426500E+00    4.4967000E-01   -2.7299000E-01
   3.5547400E+00    5.9588000E-01   -3.6771000E-01
   3.5400200E+00    7.0611000E-01    0.0000000E+00
   4.2279000E+00    5.3226000E-01    0.0000000E+00
   4.8387500E+00    2.5076000E-01    0.0000000E+00
   3.5550600E+00    3.7059000E-01   -5.9401000E-01
   2.8492200E+00    4.3026000E-01   -6.9701000E-01
   2.8491700E+00    6.9697000E-01   -4.3021000E-01
   2.1404800E+00    4.7009000E-01   -7.6877000E-01
   1.4284700E+00    4.9948000E-01   -8.1479000E-01
   1.4282600E+00    8.1543000E-01   -4.9756000E-01
   2.1403300E+00    7.6784000E-01   -4.7067000E-01
   7.1466000E-01    5.1691000E-01   -8.4100000E-01
   7.1469000E-01    8.4232000E-01   -5.1482000E-01
   7.1311000E-01    9.8871000E-01    0.0000000E+00
   1.4258100E+00    9.5800000E-01    0.0000000E+00
   2.1366800E+00    9.0283000E-01    0.0000000E+00
   2.8440800E+00    8.2123000E-01    0.0000000E+00
  -4.8267600E+00    0.0000000E+00    2.5997000E-01
  -4.9673800E+00    0.0000000E+00    1.1169000E-01
  -4.8324000E+00   -2.1791000E-01    1.3313000E-01
  -4.7182400E+00   -1.8171000E-01    2.7587000E-01
  -4.2410100E+00   -2.8370000E-01    4.4400000E-01
  -4.2278100E+00    0.0000000E+00    5.3229000E-01
  -4.2426500E+00   -4.4967000E-01    2.7299000E-01
  -3.5547400E+00   -5.9588000E-01    3.6771000E-01
  -3.5550600E+00   -3.7059000E-01    5.9401000E-01
  -2.8491700E+00   -6.9697000E-01    4.3021000E-01
  -2.8492200E+00   -4.3026000E-01    6.9701000E-01
  -3.5407600E+00    0.0000000E+00    7.0595000E-01
  -2.8438200E+00    0.0000000E+00    8.2126000E-01
  -2.1404800E+00   -4.7009000E-01    7.6877000E-01
  -2.1367200E+00    0.0000000E+00    9.0282000E-01
  -2.1403300E+00   -7.6784000E-01    4.7067000E-01
  -1.4282600E+00   -8.1543000E-01    4.9756000E-01
  -1.4284700E+00   -4.9948000E-01    8.1479000E-01
  -1.4260000E+00    0.0000000E+00    9.5798000E-01
  -7.1466000E-01   -5.1691000E-01    8.4100000E-01
  -7.1305000E-01    0.0000000E+00    9.8871000E-01
  -7.1469000E-01   -8.4232000E-01    5.1482000E-01
   0.0000000E+00   -8.5971000E-01    5.0827000E-01
   0.0000000E+00   -5.1135000E-01    8.5791000E-01
   0.0000000E+00    0.0000000E+00    1.0000000E+00
   4.7182400E+00   -1.8171000E-01    2.7587000E-01
   4.8324000E+00   -2.1791000E-01    1.3313000E-01
   4.9673800E+00    0.0000000E+00    1.1169000E-01
   4.8267600E+00    0.0000000E+00    2.5997000E-01
   4.2278100E+00    0.0000000E+00    5.3229000E-01
   4.2410100E+00   -2.8370000E-01    4.4400000E-01
   4.2426500E+00   -4.4967000E-01    2.7299000E-01
   3.5547400E+00   -5.9588000E-01    3.6771000E-01
   3.5550600E+00   -3.7059000E-01    5.9401000E-01
   2.8492200E+00   -4.3026000E-01    6.9701000E-01
   2.8491700E+00   -6.9697000E-01    4.3021000E-01
   3.5407600E+00    0.0000000E+00    7.0595000E-01
   2.8438200E+00    0.0000000E+00    8.2126000E-01
   2.1367200E+00    0.0000000E+00    9.0282000E-01
   2.1404800E+00   -4.7009000E-01    7.6877000E-01
   1.4284700E+00   -4.9948000E-01    8.1479000E-01
   1.4282600E+00   -8.1543000E-01    4.9756000E-01
   2.1403300E+00   -7.6784000E-01    4.7067000E-01
   1.4260000E+00    0.0000000E+00    9.5798000E-01
   7.1305000E-01    0.0000000E+00    9.8871000E-01
   7.1466000E-01   -5.1691000E-01    8.4100000E-01
   7.1469000E-01   -8.4232000E-01    5.1482000E-01
  -4.7182400E+00    1.8171000E-01    2.7587000E-01
  -4.8324000E+00    2.1791000E-01    1.3313000E-01
  -4.2410100E+00    2.8370000E-01    4.4400000E-01
  -4.2426500E+00    4.4967000E-01    2.7299000E-01
  -3.5547400E+00    5.9588000E-01    3.6771000E-01
  -3.5550600E+00    3.7059000E-01    5.9401000E-01
  -2.8492200E+00    4.3026000E-01    6.9701000E-01
  -2.8491700E+00    6.9697000E-01    4.3021000E-01
  -2.1404800E+00    4.7009000E-01    7.6877000E-01
  -1.4284700E+00    4.9948000E-01    8.1479000E-01
  -1.4282600E+00    8.1543000E-01    4.9756000E-01
  -2.1403300E+00    7.6784000E-01    4.7067000E-01
  -7.1466000E-01    5.1691000E-01    8.4100000E-01
   0.0000000E+00    5.1135000E-01    8.5791000E-01
   0.0000000E+00    8.5971000E-01    5.0827000E-01
  -7.1469000E-01    8.4232000E-01    5.1482000E-01
   4.8324000E+00    2.1791000E-01    1.3313000E-01
   4.7182400E+00    1.8171000E-01    2.7587000E-01
   4.2410100E+00    2.8370000E-01    4.4400000E-01
   4.2426500E+00    4.4967000E-01    2.7299000E-01
   3.5547400E+00    5.9588000E-01    3.6771000E-01
   3.5550600E+00    3.7059000E-01    5.9401000E-01
   2.8491700E+00    6.9697000E-01    4.3021000E-01
   2.8492200E+00    4.3026000E-01    6.9701000E-01
   2.1404800E+00    4.7009000E-01    7.6877000E-01
   2.1403300E+00    7.6784000E-01    4.7067000E-01
   1.4282600E+00    8.1543000E-01    4.9756000E-01
   1.4284700E+00    4.9948000E-01    8.1479000E-01
   7.1466000E-01    5.1691000E-01    8.4100000E-01
   7.1469000E-01    8.4232000E-01    5.1482000E-01
  -4.7191800E+00   -3.1896200E+00   -2.6947000E-01
  -4.8346300E+00   -3.2196400E+00   -1.2677000E-01
  -4.9677300E+00   -3.0000000E+00   -1.1107000E-01
  -4.8278700E+00   -3.0000000E+00   -2.5912000E-01
  -4.2389800E+00   -3.2884400E+00   -4.4189000E-01
  -4.2403800E+00   -3.4538800E+00   -2.6760000E-01
  -4.2258100E+00   -3.0000000E+00   -5.3299000E-01
  -3.5389000E+00   -3.0000000E+00   -7.0636000E-01
  -3.5533800E+00   -3.3733700E+00   -5.9275000E-01
  -2.8476100E+00   -3.4312300E+00   -6.9671000E-01
  -2.8476700E+00   -3.6982700E+00   -4.2858000E-01
  -3.5529700E+00   -3.5964800E+00   -3.6763000E-01
  -2.8421900E+00   -3.0000000E+00   -8.2150000E-01
  -2.1351900E+00   -3.0000000E+00   -9.0296000E-01
  -2.1393700E+00   -3.4719800E+00   -7.6755000E-01
  -1.4275400E+00   -3.4990600E+00   -8.1505000E-01
  -1.4275600E+00   -3.8152500E+00   -4.9793000E-01
  -2.1394900E+00   -3.7689400E+00   -4.6963000E-01
  -1.4252800E+00   -3.0000000E+00   -9.5803000E-01
  -7.1250000E-01   -3.0000000E+00   -9.8872000E-01
  -7.1440000E-01   -3.5136600E+00   -8.4295000E-01
   0.0000000E+00   -3.5081300E+00   -8.5979000E-01
   0.0000000E+00   -3.8585300E+00   -5.1030000E-01
  -7.1439000E-01   -3.8422400E+00   -5.1496000E-01
   0.0000000E+00   -3.0000000E+00   -1.0000000E+00
   0.0000000E+00   -4.0000000E+00    0.0000000E+00
  -7.1250000E-01   -3.9887200E+00    0.0000000E+00
  -1.4251100E+00   -3.9580400E+00    0.0000000E+00
  -2.1356600E+00   -3.9029200E+00    0.0000000E+00
  -2.8424700E+00   -3.8214600E+00    0.0000000E+00
  -3.5379900E+00   -3.7065600E+00    0.0000000E+00
  -4.2249800E+00   -3.5332400E+00    0.0000000E+00
  -4.8411700E+00   -3.2489000E+00    0.0000000E+00
  -5.0000000E+00   -3.0000000E+00    0.0000000E+00
   4.8278700E+00   -3.0000000E+00   -2.5912000E-01
   4.9677300E+00   -3.0000000E+00   -1.1107000E-01
   4.8346300E+00   -3.2196400E+00   -1.2677000E-01
   4.7191800E+00   -3.1896200E+00   -2.6947000E-01
   4.2403800E+00   -3.4538800E+00   -2.6760000E-01
   4.2389800E+00   -3.2884400E+00   -4.4189000E-01
   4.2258100E+00   -3.0000000E+00   -5.3299000E-01
   3.5533800E+00   -3.3733700E+00   -5.9275000E-01
   3.5389000E+00   -3.0000000E+00   -7.0636000E-01
   3.5529700E+00   -3.5964800E+00   -3.6763000E-01
   2.8476700E+00   -3.6982700E+00   -4.2858000E-01
   2.8476100E+00   -3.4312300E+00   -6.9671000E-01
   2.8421900E+00   -3.0000000E+00   -8.2150000E-01
   2.1393700E+00   -3.4719800E+00   -7.6755000E-01
   2.1351900E+00   -3.0000000E+00   -9.0296000E-01
   2.1394900E+00   -3.7689400E+00   -4.6963000E-01
   1.4275600E+00   -3.8152500E+00   -4.9793000E-01
   1.4275400E+00   -3.4990600E+00   -8.1505000E-01
   1.4252800E+00   -3.0000000E+00   -9.5803000E-01
   7.1440000E-01   -3.5136600E+00   -8.4295000E-01
   7.1250000E-01   -3.0000000E+00   -9.8872000E-01
   7.1439000E-01   -3.8422400E+00   -5.1496000E-01
   7.1250000E-01   -3.9887200E+00    0.0000000E+00
   1.4251100E+00   -3.9580400E+00    0.0000000E+00
   2.1356600E+00   -3.9029200E+00    0.0000000E+00
   2.8424700E+00   -3.8214600E+00    0.0000000E+00
   3.5379900E+00   -3.7065600E+00    0.0000000E+00
   4.2249800E+00   -3.5332400E+00    0.0000000E+00
   4.8411700E+00   -3.2489000E+00    0.0000000E+00
   5.0000000E+00   -3.0000000E+00    0.0000000E+00
  -4.8346300E+00   -2.7803600E+00   -1.2677000E-01
  -4.7191800E+00   -2.8103800E+00   -2.6947000E-01
  -4.2403800E+00   -2.5461200E+00   -2.6760000E-01
  -4.2389800E+00   -2.7115600E+00   -4.4189000E-01
  -3.5533800E+00   -2.6266300E+00   -5.9275000E-01
  -3.5529700E+00   -2.4035200E+00   -3.6763000E-01
  -2.8476700E+00   -2.3017300E+00   -4.2858000E-01
  -2.8476100E+00   -2.5687700E+00   -6.9671000E-01
  -2.1393700E+00   -2.5280200E+00   -7.6755000E-01
  -2.1394900E+00   -2.2310600E+00   -4.6963000E-01
  -1.4275600E+00   -2.1847500E+00   -4.9793000E-01
  -1.4275400E+00   -2.5009400E+00   -8.1505000E-01
  -7.1440000E-01   -2.4863400E+00   -8.4295000E-01
  -7.1439000E-01   -2.1577600E+00   -5.1496000E-01
   0.0000000E+00   -2.1414700E+00   -5.1030000E-01
   0.0000000E+00   -2.4918700E+00   -8.5979000E-01
  -7.1250000E-01   -2.0112800E+00    0.0000000E+00
   0.0000000E+00   -2.0000000E+00    0.0000000E+00
  -1.4251100E+00   -2.0419600E+00    0.0000000E+00
  -2.1356600E+00   -2.0970800E+00    0.0000000E+00
  -2.8424700E+00   -2.1785400E+00    0.0000000E+00
  -3.5379900E+00   -2.2934400E+00    0.0000000E+00
  -4.2249800E+00   -2.4667600E+00    0.0000000E+00
  -4.8411700E+00   -2.7511000E+00    0.0000000E+00
   4.7191800E+00   -2.8103800E+00   -2.6947000E-01
   4.8346300E+00   -2.7803600E+00   -1.2677000E-01
   4.2389800E+00   -2.7115600E+00   -4.4189000E-01
   4.2403800E+00   -2.5461200E+00   -2.6760000E-01
   3.5533800E+00   -2.6266300E+00   -5.9275000E-01
   2.8476100E+00   -2.5687700E+00   -6.9671000E-01
   2.8476700E+00   -2.3017300E+00   -4.2858000E-01
   3.5529700E+00   -2.4035200E+00   -3.6763000E-01
   2.1393700E+00   -2.5280200E+00   -7.6755000E-01
   1.4275400E+00   -2.5009400E+00   -8.1505000E-01
   1.4275600E+00   -2.1847500E+00   -4.9793000E-01
   2.1394900E+00   -2.2310600E+00   -4.6963000E-01
   7.1440000E-01   -2.4863400E+00   -8.4295000E-01
   7.1439000E-01   -2.1577600E+00   -5.1496000E-01
   7.1250000E-01   -2.0112800E+00    0.0000000E+00
   1.4251100E+00   -2.0419600E+00    0.0000000E+00
   2.1356600E+00   -2.0970800E+00    0.0000000E+00
   2.8424700E+00   -2.1785400E+00    0.0000000E+00
   3.5379900E+00   -2.2934400E+00    0.0000000E+00
   4.2249800E+00   -2.4667600E+00    0.0000000E+00
   4.8411700E+00   -2.7511000E+00    0.0000000E+00
  -4.8278700E+00   -3.0000000E+00    2.5912000E-01
  -4.9677300E+00   -3.0000000E+00    1.1107000E-01
  -4.8346300E+00   -3.2196400E+00    1.2677000E-01
  -4.7191800E+00   -3.1896200E+00    2.6947000E-01
  -4.2403800E+00   -3.4538800E+00    2.6760000E-01
  -4.2389800E+00   -3.2884400E+00    4.4189000E-01
  -4.2258100E+00   -3.0000000E+00    5.3299000E-01
  -3.5533800E+00   -3.3733700E+00    5.9275000E-01
  -3.5389000E+00   -3.0000000E+00    7.0636000E-01
  -3.5529700E+00   -3.5964800E+00    3.6763000E-01
  -2.8476700E+00   -3.6982700E+00    4.2858000E-01
  -2.8476100E+00   -3.4312300E+00    6.9671000E-01
  -2.8421900E+00   -3.0000000E+00    8.2150000E-01
  -2.1393700E+00   -3.4719800E+00    7.6755000E-01
  -2.1351900E+00   -3.0000000E+00    9.0296000E-01
  -2.1394900E+00   -3.7689400E+00    4.6963000E-01
  -1.4275600E+00   -3.8152500E+00    4.9793000E-01
  -1.4275400E+00   -3.4990600E+00    8.1505000E-01
  -1.4252800E+00   -3.0000000E+00    9.5803000E-01
  -7.1440000E-01   -3.5136600E+00    8.4295000E-01
  -7.1250000E-01   -3.0000000E+00    9.8872000E-01
  -7.1439000E-01   -3.8422400E+00    5.1496000E-01
   0.0000000E+00   -3.8585300E+00    5.1030000E-01
   0.0000000E+00   -3.5081300E+00    8.5979000E-01
   0.0000000E+00   -3.0000000E+00    1.0000000E+00
   4.7191800E+00   -3.1896200E+00    2.6947000E-01
   4.8346300E+00   -3.2196400E+00    1.2677000E-01
   4.9677300E+00   -3.0000000E+00    1.1107000E-01
   4.8278700E+00   -3.0000000E+00    2.5912000E-01
   4.2389800E+00   -3.2884400E+00    4.4189000E-01
   4.2403800E+00   -3.4538800E+00    2.6760000E-01
   4.2258100E+00   -3.0000000E+00    5.3299000E-01
   3.5389000E+00   -3.0000000E+00    7.0636000E-01
   3.5533800E+00   -3.3733700E+00    5.9275000E-01
   2.8476100E+00   -3.4312300E+00    6.9671000E-01
   2.8476700E+00   -3.6982700E+00    4.2858000E-01
   3.5529700E+00   -3.5964800E+00    3.6763000E-01
   2.8421900E+00   -3.0000000E+00    8.2150000E-01
   2.1351900E+00   -3.0000000E+00    9.0296000E-01
   2.1393700E+00   -3.4719800E+00    7.6755000E-01
   1.4275400E+00   -3.4990600E+00    8.1505000E-01
   1.4275600E+00   -3.8152500E+00    4.9793000E-01
   2.1394900E+00   -3.7689400E+00    4.6963000E-01
   1.4252800E+00   -3.0000000E+00    9.5803000E-01
   7.1250000E-01   -3.0000000E+00    9.8872000E-01
   7.1440000E-01   -3.5136600E+00    8.4295000E-01
   7.1439000E-01   -3.8422400E+00    5.1496000E-01
  -4.7191800E+00   -2.8103800E+00    2.6947000E-01
  -4.8346300E+00   -2.7803600E+00    1.2677000E-01
  -4.2389800E+00   -2.7115600E+00    4.4189000E-01
  -4.2403800E+00   -2.5461200E+00    2.6760000E-01
  -3.5533800E+00   -2.6266300E+00    5.9275000E-01
  -2.8476100E+00   -2.5687700E+00    6.9671000E-01
  -2.8476700E+00   -2.3017300E+00    4.2858000E-01
  -3.5529700E+00   -2.4035200E+00    3.6763000E-01
  -2.1393700E+00   -2.5280200E+00    7.6755000E-01
  -1.4275400E+00   -2.5009400E+00    8.1505000E-01
  -1.4275600E+00   -2.1847500E+00    4.9793000E-01
  -2.1394900E+00   -2.2310600E+00    4.6963000E-01
  -7.1440000E-01   -2.4863400E+00    8.4295000E-01
   0.0000000E+00   -2.4918700E+00    8.5979000E-01
   0.0000000E+00   -2.1414700E+00    5.1030000E-01
  -7.1439000E-01   -2.1577600E+00    5.1496000E-01
   4.8346300E+00   -2.7803600E+00    1.2677000E-01
   4.7191800E+00   -2.8103800E+00    2.6947000E-01
   4.2403800E+00   -2.5461200E+00    2.6760000E-01
   4.2389800E+00   -2.7115600E+00    4.4189000E-01
   3.5533800E+00   -2.6266300E+00    5.9275000E-01
   3.5529700E+00   -2.4035200E+00    3.6763000E-01
   2.8476700E+00   -2.3017300E+00    4.2858000E-01
   2.8476100E+00   -2.5687700E+00    6.9671000E-01
   2.1393700E+00   -2.5280200E+00    7.6755000E-01
   2.1394900E+00   -2.2310600E+00    4.6963000E-01
   1.4275600E+00   -2.1847500E+00    4.9793000E-01
   1.4275400E+00   -2.5009400E+00    8.1505000E-01
   7.1440000E-01   -2.4863400E+00    8.4295000E-01
   7.1439000E-01   -2.1577600E+00    5.1496000E-01
    1          2          3          4
    1          4          5          6
    1          6          7          2
    7          8          9         10
    7         10         11          2
    7          6         12          8
   12         13         14          8
   12          6          5         15
   12         15         16         13
   16         17         18         13
   18         19         20         21
   18         21         14         13
   18         17         22         19
   22         23         24         19
   24         25         26         27
   24         27         20         19
   24         23         28         25
   27         26         29         30
   27         30         31         20
   31         32         21         20
   32         33         14         21
   33          9          8         14
   11         34          3          2
   35         36         37         38
   39         40         35         38
   37         41         39         38
   42         43         44         41
   37         45         42         41
   44         46         39         41
   44         47         48         46
   49         40         39         46
   48         50         49         46
   48         51         52         50
   53         54         55         51
   48         47         53         51
   55         56         52         51
   55         57         58         56
   59         26         25         57
   55         54         59         57
   25         28         58         57
   60         29         26         59
   54         61         60         59
   54         53         62         61
   53         47         63         62
   47         44         43         63
   37         36         64         45
    4          3         65         66
   67          5          4         66
   65         68         67         66
   69         70         71         68
   65         72         69         68
   71         73         67         68
   71         74         75         73
   15          5         67         73
   75         16         15         73
   75         76         17         16
   77         78         79         76
   75         74         77         76
   79         22         17         76
   79         80         23         22
   81         82         83         80
   79         78         81         80
   83         28         23         80
   84         85         82         81
   78         86         84         81
   78         77         87         86
   77         74         88         87
   74         71         70         88
   65          3         34         72
   89         90         36         35
   89         35         40         91
   89         91         92         90
   92         93         94         95
   92         95         96         90
   92         91         97         93
   97         98         99         93
   97         91         40         49
   97         49         50         98
   50         52        100         98
  100        101        102        103
  100        103         99         98
  100         52         56        101
   56         58        104        101
  104         83         82        105
  104        105        102        101
  104         58         28         83
  105         82         85        106
  105        106        107        102
  107        108        103        102
  108        109         99        103
  109         94         93         99
   96         64         36         90
  110        111        112        113
  114        115        110        113
  112        116        114        113
   10          9        117        116
  112         11         10        116
  117        118        114        116
  117        119        120        118
  121        115        114        118
  120        122        121        118
  120        123        124        122
  125        126        127        123
  120        119        125        123
  127        128        124        123
  127        129        130        128
  131        132        133        129
  127        126        131        129
  133        134        130        129
   30         29        132        131
  126         31         30        131
  126        125         32         31
  125        119         33         32
  119        117          9         33
  112        111         34         11
  135        136        137        138
  135        138        139        140
  135        140        141        136
  141        142         43         42
  141         42         45        136
  141        140        143        142
  143        144        145        142
  143        140        139        146
  143        146        147        144
  147        148        149        144
  149        150        151        152
  149        152        145        144
  149        148        153        150
  153        154        155        150
  155        133        132        156
  155        156        151        150
  155        154        134        133
  156        132         29         60
  156         60         61        151
   61         62        152        151
   62         63        145        152
   63         43        142        145
   45         64        137        136
  157        158        111        110
  157        110        115        159
  157        159        160        158
  160        161         70         69
  160         69         72        158
  160        159        162        161
  162        163        164        161
  162        159        115        121
  162        121        122        163
  122        124        165        163
  165        166        167        168
  165        168        164        163
  165        124        128        166
  128        130        169        166
  169        170        171        172
  169        172        167        166
  169        130        134        170
  172        171         85         84
  172         84         86        167
   86         87        168        167
   87         88        164        168
   88         70        161        164
   72         34        111        158
  138        137        173        174
  175        139        138        174
  173        176        175        174
   95         94        177        176
  173         96         95        176
  177        178        175        176
  177        179        180        178
  146        139        175        178
  180        147        146        178
  180        181        148        147
  182        183        184        181
  180        179        182        181
  184        153        148        181
  184        185        154        153
  186        171        170        185
  184        183        186        185
  170        134        154        185
  106         85        171        186
  183        107        106        186
  183        182        108        107
  182        179        109        108
  179        177         94        109
  173        137         64         96
  187        188        189        190
  187        191        192        188
  187        190        193        191
  193        194        195        191
  195        196        197        198
  195        198        192        191
  195        194        199        196
  199        200        201        196
  201        202        203        204
  201        204        197        196
  201        200        205        202
  205        206        207        202
  207        208        209        210
  207        210        203        202
  207        206        211        208
  210        209        212        213
  210        213        214        203
  214        215        204        203
  215        216        197        204
  216        217        198        197
  217        218        192        198
  218        219        188        192
  219        220        189        188
  221        222        223        224
  223        225        226        224
  226        227        221        224
  226        228        229        227
  230        231        232        228
  226        225        230        228
  232        233        229        228
  232        234        235        233
  236        237        238        234
  232        231        236        234
  238        239        235        234
  238        240        241        239
  242        209        208        240
  238        237        242        240
  208        211        241        240
  243        212        209        242
  237        244        243        242
  237        236        245        244
  236        231        246        245
  231        230        247        246
  230        225        248        247
  225        223        249        248
  223        222        250        249
  190        189        251        252
  251        253        254        252
  254        193        190        252
  254        255        194        193
  256        257        258        255
  254        253        256        255
  258        199        194        255
  258        259        200        199
  260        261        262        259
  258        257        260        259
  262        205        200        259
  262        263        206        205
  264        265        266        263
  262        261        264        263
  266        211        206        263
  267        268        265        264
  261        269        267        264
  261        260        270        269
  260        257        271        270
  257        256        272        271
  256        253        273        272
  253        251        274        273
  251        189        220        274
  275        276        222        221
  275        277        278        276
  275        221        227        277
  227        229        279        277
  279        280        281        282
  279        282        278        277
  279        229        233        280
  233        235        283        280
  283        284        285        286
  283        286        281        280
  283        235        239        284
  239        241        287        284
  287        266        265        288
  287        288        285        284
  287        241        211        266
  288        265        268        289
  288        289        290        285
  290        291        286        285
  291        292        281        286
  292        293        282        281
  293        294        278        282
  294        295        276        278
  295        250        222        276
  296        297        298        299
  298        300        301        299
  301        302        296        299
  301        303        304        302
  305        306        307        303
  301        300        305        303
  307        308        304        303
  307        309        310        308
  311        312        313        309
  307        306        311        309
  313        314        310        309
  313        315        316        314
  317        318        319        315
  313        312        317        315
  319        320        316        315
  213        212        318        317
  312        214        213        317
  312        311        215        214
  311        306        216        215
  306        305        217        216
  305        300        218        217
  300        298        219        218
  298        297        220        219
  321        322        323        324
  321        325        326        322
  321        324        327        325
  327        328        329        325
  329        330        331        332
  329        332        326        325
  329        328        333        330
  333        334        335        330
  335        336        337        338
  335        338        331        330
  335        334        339        336
  339        340        341        336
  341        319        318        342
  341        342        337        336
  341        340        320        319
  342        318        212        243
  342        243        244        337
  244        245        338        337
  245        246        331        338
  246        247        332        331
  247        248        326        332
  248        249        322        326
  249        250        323        322
  343        344        297        296
  343        345        346        344
  343        296        302        345
  302        304        347        345
  347        348        349        350
  347        350        346        345
  347        304        308        348
  308        310        351        348
  351        352        353        354
  351        354        349        348
  351        310        314        352
  314        316        355        352
  355        356        357        358
  355        358        353        352
  355        316        320        356
  358        357        268        267
  358        267        269        353
  269        270        354        353
  270        271        349        354
  271        272        350        349
  272        273        346        350
  273        274        344        346
  274        220        297        344
  324        323        359        360
  359        361        362        360
  362        327        324        360
  362        363        328        327
  364        365        366        363
  362        361        364        363
  366        333        328        363
  366        367        334        333
  368        369        370        367
  366        365        368        367
  370        339        334        367
  370        371        340        339
  372        357        356        371
  370        369        372        371
  356        320        340        371
  289        268        357        372
  369        290        289        372
  369        368        291        290
  368        365        292        291
  365        364        293        292
  364        361        294        293
  361        359        295        294
  359        323        250        295

~~都看到这了，点个赞再走吧~~~