Paste from 11 at 2022-07-14 10:00:19I want to paste

program main1
 implicit  none
 real (kind=4) :: r,s,k,p,pc,pz
 real (kind=8) :: u,v,q,x
 integer , parameter :: L=5 
 integer :: t
 integer :: status=0
 integer :: pipe
 integer :: pip !用来赋值管道序号
 integer :: time
 real , allocatable :: qo(:,:)
 real , allocatable :: q1(:,:)
 real , allocatable :: xl(:,:)
 real , allocatable :: u1(:,:)
 real , allocatable :: v1(:,:)
 real , allocatable :: qw(:,:)
 integer , allocatable :: timemax(:)
 real , allocatable :: sl(:)
 real (kind=4) :: ansq
 real (kind=4) :: ansqo
 real (kind=4) :: ansqw
 real (kind=4) :: ww
 real (kind=4) :: ko
 real (kind=4) :: kw 
 real (kind=4) :: ansx
 real (kind=4) :: sw
 real (kind=4) :: vw
 real (kind=4) :: vz
 real (kind=4) :: ansvw
 real (kind=4) :: ansvz
logical alive
 write(*,*) "压差p(MPa)"
 read(*,*) p
 inquire(file="3.txt",exist=alive)
 
 if(.not.alive)then
     write(*,*) "文件不存在"
     stop
 end if
 
 open(unit=10,file="3.txt")
 pipe=0
 
 do while (.true.)
     read(10,"(F9.3)",iostat=status)r
     if(status/=0) exit
     pipe=pipe+1
 end do
 
 rewind(10)
 !此处读取管道个数的程序，即为pipe的大小
    write(*,*)"pipe=" , pipe
    allocate(timemax(pipe))
    allocate(sl(pipe))
 do pip = 1 , pipe
  
   read(10,*,iostat=status)r
   pc=2.0*30.0*0.000001/(r*0.001)
   pz=p-pc
   if  (pz>0) then 
       write(*,*) "pz=" , pz 
   s=3.14*(r*0.0001/2.0)**2
   x=0
   t=0
     do while (x<=L)
       u=(5.0*(L-x)**2+0.75*x**2)/((L-x)**2+x**2)
       v=((r*0.000001/2.0))**2*pz*1000000/(8.0*u*0.001*0.05*100.0)
       q=v*s
       x=x+v
       t=t+1
       timemax(pip)=t
    write(*,*)"x=" , x
    write(*,*)"pip=" , pip
   
     end do
   else
       write(*,*) "无法驱替"
       GOTO 100
   end if  
 end do
100     write(*,*)"t=" , t
     timemax(pipe)=t
     write(*,*)"timemax=" , t
     time=t
     
     allocate(qo(pipe,time))
     allocate(q1(pipe,time))
     allocate(qw(pipe,time))
     allocate(xl(pipe,time))
     allocate(u1(pipe,time))
     allocate(v1(pipe,time))
     
 !计算某一时刻含水率
   
    xl(:,:)=0
   
    do t = 1 , timemax(pipe)
        
        rewind(10)
   do pip = 1 , pipe
       
      
       read(10,*,iostat=status)r
       sl(pip)=3.14*(r*0.0001/2.0)**2
       pc=2.0*30.0*0.000001/(r*0.001)
       pz=p-pc
       
        if  (pz>0) then
             
            
                 ansq=0
                 ansqo=0
                 ansqw=0
                 ansx=0
             if (t<timemax(pip)) then
                     write(*,*) "pip=" , pip
                     write(*,*) "t=" , t
                     
                
                 
                     u1(pip,t)=(5.0*(L-xl(pip,t))**2+0.75*xl(pip,t)**2)/((L-xl(pip,t))**2+xl(pip,t)**2)
                     v1(pip,t)=((r*0.000001/2.0))**2*pz*1000000/(8*u1(pip,t)*0.001*0.05*100.0)
                     q1(pip,t)=v1(pip,t)*sl(pip)
                     xl(pip,t+1)=xl(pip,t)+v1(pip,t)
                     qw(pip,t)=0
                     qo(pip,t)=q1(pip,t)
                     
                 
             else
                 xl(pip,t)= L
                 u1(pip,t)=0.75
                 v1(pip,t)=((r*0.000001/2.0))**2*p*1000000/(8.0*u1(pip,t)*0.001*0.05*100.0)
                 q1(pip,t)=v1(pip,t)*sl(pip)
                 qw(pip,t)=q1(pip,t)
                 qo(pip,t)=0

             end if
             
           write(*,*) "t=",t
           write(*,*) "pip=",pip
           write(*,*) "u1(pip,t)=" , u1(pip,t)
           write(*,*) "qo(pip,t)=" , qo(pip,t)
           write(*,*) "qw(pip,t)=" , qw(pip,t)
        else
           
        write(*,*) "无法驱替"
        xl(pip,t)= 0
        u1(pip,t)=5
        v1(pip,t)=0
        q1(pip,t)=0
        qw(pip,t)=0
        qo(pip,t)=0
    
       end if 
  end do
       do pip = 1 , pipe
  
          ansq=ansq+qo(pip,t)+qw(pip,t)
          ansqo=qo(pip,t)+ansqo
          ansqw=qw(pip,t)+ansqw
          ww=ansqw/ansq
         
end do
  write(*,*) "ww=",ww

    end do
    t=t-1
    ansvw=0
    ansvz=0
  do pip = 1 , pipe
    write(*,*) "t=" , t
    write(*,*) "pip=" , pip
   
    vw=xl(pip,t)*sl(pip)
    ansvw=ansvw+vw
    vz=5*sl(pip)
    ansvz=ansvz+vz
    write(*,*) "vw=" , vw
    write(*,*) "vz=" , vz
  end do
    write(*,*) "ansvw=" , ansvw
    write(*,*) "ansvz=" , ansvz
    sw=ansvw/ansvz
    write(*,*) "sw=",sw
if (ww>0.98) GOTO 200
200   write(*,*) "t=",t
   write(*,*) "ansq=",ansq
   write(*,*) "ansqo=",ansqo
   write(*,*) "ansqw=",ansqw
   !计算相渗透率
   
   ko=(ansqo*5.0*L)/(1.5625*p)
   kw=(ansqw*0.75*L)/(1.5625*p)
   open(unit=11,file="1.txt",position="append")
   write(11,*) "ww=",ww
   write(11,*) "ko=",ko
   write(11,*) "kw=",kw
   write(11,*) "sw=",sw
   write(*,*) "ww=",ww
   write(*,*) "ko=",ko
   write(*,*) "kw=",kw
   write(*,*) "sw=",sw
   close(unit=11)    
   
 
   
       

 
STOP
end