SimulationCore2/FluidFlow/Well_Pressure_Data_Transfer.f90

SUBROUTINE WellPressureDataTransfer
    
    !! This subroutine calculates pressure distribution in well
    !! We divide well hydraulic system to 4 zones. 1- pump to string zone; this zone is horizontal
    !! 2- string zone; this zone includes drill pipes and collars and is vertical or neraly vertical
    !! 3- annulus zone; this zone includes annulus space between string and openhole or casing
    !! 4- Openhole zone; this zone includes spaces under bit to bottom of the wellbore
    !! The last two zone includes annulus space and open hole, are influenced by kicks and have crucial role in 
    !! computing pressure distribution because reference point (back pressure in well head open condition 
    !! and bottomhole pressure in well head closed condition) are in one of them. because of this, pressure
    !! distribution calculates in two subroutines : PressureAnnAndOHDistribution and PressurePumptoBitDistribution
    !! pressure change through well consist of pressure change due to friction, change of elevation and acceleration
    !! frictional pressure drop calculated by subroutine Sub_Press_Drop_Calculator
    !! This subroutine includes hydrostatic pressure and later will include effect of corss-section area change
    !! which causes acceleration or decceleration of flow
    
    USE FricPressDropVarsModule
    use MudSystemModule
    USE MudSystemVARIABLES
use SimulationVariables !@@@
    USE GeoElements_FluidModule
    USE Fluid_Flow_Startup_Vars
    use KickVARIABLESModule
    USE CMudPropertiesVariables
    USE CBopStackVariables
    USE sROP_Variables
    USE CHOKEVARIABLES
use SimulationVariables !@, PI2 => PI
    use CChokeManifoldVariables
    use SimulationVariables
    USE CBitProblemsVariables
    USE CChokeProblemsVariables
    USE CStringConfigurationVariables
    use SimulationVariables !@
    use CDataDisplayConsole
    use SimulationVariables
    use SimulationVariables
    
    IMPLICIT NONE
    
    INTEGER :: i , j
    REAL :: InstantaneousTotalOpenChokeArea , VolumeDensityProduct
    
    KickVARIABLES%ChokeIsClosing = .FALSE.
    data%State%FricPressDrop%Kchoke = 0.0
    data%State%FricPressDrop%KBit = 0.0
    
    
    !!!!!!!!!!!!!!!! bit is present
    BitTrue = data%Configuration%StringConfiguration%StringConfigurations(1)%ComponentType  ==  Bit_ComponentType           ! Bit_ComponentType = 0
    BitCd = 0.98
    BitNozzDia = 32.0 * data%Configuration%StringConfiguration%BitDefinition%BitNozzleSize       ! nozzle diameter in 1/32 in
    BitNozzleArea = 7.6699E-4 * BitNozzDia**2
    BitNozzleNum = data%Configuration%StringConfiguration%BitDefinition%BitNozzleNo - data%State%FricPressDrop%BitJetsPlugged * data%problems%BitProblems%PlugJetsCount
    BitTotNozzArea = BitNozzleNum * BitNozzleArea + 0.5 * data%State%FricPressDrop%BitJetsWashedOut * data%problems%BitProblems%JetWashoutCount * BitNozzleArea
    IF (BitNozzleNum == 0) THEN
        data%State%FricPressDrop%BitTotallyPluged = .TRUE.
    ELSE
        data%State%FricPressDrop%BitTotallyPluged = .FALSE.
    END IF
    
    data%State%FricPressDrop%FloatValveIn = data%Configuration%StringConfiguration%BitDefinition%FloatValve
    !FloatValveOpen = .TRUE.
    !IF (JetWashoutCount > 0 .OR. PlugJetsCount > 0) THEN
    !    WRITE (*,*) 'PlugJetsCount' , PlugJetsCount
    !    WRITE (*,*) 'JetWashoutCount', JetWashoutCount
    !END IF
    
    !WRITE (*,*) 'BitNozzleNum', BitNozzleNum
    !WRITE (*,*) 'BitTotNozzArea', BitTotNozzArea
    !WRITE (*,*) 'Float Valve??', BitDefinition%FloatValve
    
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    IF(ALLOCATED(KickVARIABLES%GasPocketFlowEl))      DEALLOCATE(KickVARIABLES%GasPocketFlowEl)
    
!!!!!!!!!!!!!!!!!!!!  Well Head Condition ( Open or Closed )
    KickVARIABLES%ChokeKroneckerDelta = 0

    
    IF (data%State%MudSystem%WellToPitsOpen == .FALSE. .AND. data%State%MudSystem%WellToChokeManifoldOpen == .FALSE.) THEN
        
        KickVARIABLES%WellHeadOpen = .FALSE.
        !WRITE (*,*) ' Well Head is closed '
    ELSE 
        KickVARIABLES%WellHeadOpen = .TRUE.
    END IF
    
    IF (data%State%MudSystem%WellToPitsOpen == .FALSE. .AND. data%State%MudSystem%WellToChokeManifoldOpen)        KickVARIABLES%ChokeKroneckerDelta = 1
    !WRITE (*,*) ' WelltoPitsOpen=' , WelltoPitsOpen, 'WellToChokeManifoldOpen= ' , WellToChokeManifoldOpen
    
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    
    data%State%FricPressDrop%NoHorizontalEl = data%State%MudSystem%NoHorizontalMudElements
    data%State%FricPressDrop%NoStringEl = data%State%MudSystem%NoStringMudElements
    data%State%FricPressDrop%NoAnnulusEl = data%State%MudSystem%NoCasingMudElements
    data%State%FricPressDrop%NoWellToChokeEl = 0
    IF (KickVARIABLES%ChokeKroneckerDelta == 1)       data%State%FricPressDrop%NoWellToChokeEl = data%State%MudSystem%ChokeLine_Density%Length()
    data%State%FricPressDrop%NoOpenHoleEl = data%State%MudSystem%NoBottomHoleMudElements
    data%State%FricPressDrop%NumbEl = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + data%State%FricPressDrop%NoAnnulusEl + data%State%FricPressDrop%NoWellToChokeEl + data%State%FricPressDrop%NoOpenHoleEl
    IF(ALLOCATED(FlowEl))       DEALLOCATE(FlowEl)
    ALLOCATE(FlowEl(data%State%FricPressDrop%NumbEl))
    
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  Transfering Data from module:MudSystemVARIABLES  to  module:FricPressDropVars to calculate Pressure distribution in flow path
    !!!!!!!!!!!!!!!!!!!!!!! Horizontal line
    !WRITE (*,*) 'H', NoHorizontalEl
    DO i = 1 , data%State%FricPressDrop%NoHorizontalEl
        FlowEl(i)%StartX = data%State%MudSystem%Xstart_MudElement%Array(i) 
        FlowEl(i)%EndX = data%State%MudSystem%Xend_MudElement%Array(i)
        FlowEl(i)%StartTVD = 0.d0
        FlowEl(i)%EndTVD = 0.d0
        FlowEl(i)%Id = 0.0
        FlowEl(i)%Od = data%State%MudSystem%PipeOD_MudElement%Array(i)
        FlowEl(i)%density = data%State%MudSystem%Density_MudElement%Array(i)
        
        FlowEl(i)%Dhyd = FlowEl(i)%Od
        FlowEl(i)%Length = ABS(FlowEl(i)%EndX - FlowEl(i)%StartX)
        FlowEl(i)%DepthDiff = 0.d0
        FlowEl(i)%Volume = PI / 4.0 * FlowEl(i)%Od**2 * REAL(FlowEl(i)%Length) / Convfttoinch**2
        IF (data%State%MudSystem%MudType_MudElement%Array(i) == 0) THEN           ! = 0 for mud, = 1 for gas kick, = 4 for air
            FlowEl(i)%MaterialType = 0
        ELSE IF (data%State%MudSystem%MudType_MudElement%Array(i) == 4) THEN
            FlowEl(i)%MaterialType = 4
        ELSE
            FlowEl(i)%MaterialType = 1
        END IF
        
        FlowEl(i)%dPdLgrav = 0.0
        FlowEl(i)%StaticPressDiff = 0.0


        
    END DO
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    !!!!!!!!!!!!!!!!!!!!!!!!! String
    data%State%FricPressDrop%StringFirstEl = data%State%FricPressDrop%NoHorizontalEl + 1
    data%State%FricPressDrop%StringLastEl = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl
    DO i = data%State%FricPressDrop%NoHorizontalEl + 1 , data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl
        !WRITE (*,*) 'ST , i' , i
        FlowEl(i)%StartX = data%State%MudSystem%Xstart_MudElement%Array(i)
        !WRITE (*,*) 'StartX', FlowEl(i)%StartX
        FlowEl(i)%EndX = data%State%MudSystem%Xend_MudElement%Array(i)
        !WRITE (*,*) 'EndX', FlowEl(i)%EndX
        FlowEl(i)%StartTVD = data%State%MudSystem%TVDstart_MudElement%Array(i)
        !WRITE (*,*) 'StartTVD', FlowEl(i)%StartTVD
        FlowEl(i)%EndTVD = data%State%MudSystem%TVDend_MudElement%Array(i)
        !WRITE (*,*) 'EndTVD', FlowEl(i)%EndTVD
        FlowEl(i)%Id = 0.0
        !WRITE (*,*) 'Id', FlowEl(i)%Id
        FlowEl(i)%Od = data%State%MudSystem%PipeOD_MudElement%Array(i)
        !WRITE (*,*) 'Od', FlowEl(i)%Od
        FlowEl(i)%Density = data%State%MudSystem%Density_MudElement%Array(i)
        !WRITE (*,*) 'Density', FlowEl(i)%Density
        
        
        FlowEl(i)%Dhyd = FlowEl(i)%Od
        FlowEl(i)%Length = ABS(FlowEl(i)%EndX - FlowEl(i)%StartX)
        !WRITE (*,*) 'Length', FlowEl(i)%Length
        FlowEl(i)%DepthDiff = ABS(FlowEl(i)%StartTVD - FlowEl(i)%EndTVD)
        !WRITE (*,*) 'DepthDiff', FlowEl(i)%DepthDiff
        FlowEl(i)%Area = PI / 4.0 * FlowEl(i)%Od**2 / Convfttoinch**2
        FlowEl(i)%Volume = FlowEl(i)%Area * REAL(FlowEl(i)%Length)
        IF (data%State%MudSystem%MudType_MudElement%Array(i) == 0) THEN           ! = 0 for mud, = 1 for gas kick, = 4 for air
            FlowEl(i)%MaterialType = 0
        ELSE IF (data%State%MudSystem%MudType_MudElement%Array(i) == 4) THEN
            FlowEl(i)%MaterialType = 4
        ELSE
            FlowEl(i)%MaterialType = 1
        END IF
        
        IF (FlowEl(i)%MaterialType == 1 .OR. FlowEl(i)%MaterialType == 4) THEN
            FlowEl(i)%dPdLgrav = 0.0
        ELSE
            FlowEl(i)%dPdLgrav = 0.052 * FlowEl(i)%Density
        END IF
        FlowEl(i)%StaticPressDiff = FlowEl(i)%dPdLgrav * REAL(FlowEl(i)%DepthDiff)
        


        
    END DO
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    !!!!!!!!!!!!!!!!!!!!!!!!!!!! Annulus
    data%State%FricPressDrop%AnnulusFirstEl = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + 1
    data%State%FricPressDrop%AnnulusLastEl = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + data%State%FricPressDrop%NoAnnulusEl
    DO i = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + 1 , data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + data%State%FricPressDrop%NoAnnulusEl
        FlowEl(i)%StartX = data%State%MudSystem%Xstart_MudElement%Array(i)
        FlowEl(i)%EndX = data%State%MudSystem%Xend_MudElement%Array(i)
        FlowEl(i)%StartTVD = data%State%MudSystem%TVDstart_MudElement%Array(i)
        FlowEl(i)%EndTVD = data%State%MudSystem%TVDend_MudElement%Array(i)
        FlowEl(i)%Id = data%State%MudSystem%PipeID_MudElement%Array(i)
        FlowEl(i)%Od = data%State%MudSystem%PipeOD_MudElement%Array(i)
        FlowEl(i)%Density = data%State%MudSystem%Density_MudElement%Array(i)        
        
        FlowEl(i)%Dhyd = FlowEl(i)%Od - FlowEl(i)%Id
        FlowEl(i)%Length = ABS(FlowEl(i)%EndX - FlowEl(i)%StartX)
        FlowEl(i)%DepthDiff = ABS(FlowEl(i)%StartTVD - FlowEl(i)%EndTVD)
        FlowEl(i)%Area = PI / 4.0 * (FlowEl(i)%Od**2 - FlowEl(i)%Id**2) / Convfttoinch**2
        FlowEl(i)%Volume = FlowEl(i)%Area * REAL(FlowEl(i)%Length)
        IF (data%State%MudSystem%MudType_MudElement%Array(i) == 0) THEN           ! = 0 for mud, = 1 for gas kick, = 4 for air
            FlowEl(i)%MaterialType = 0
        ELSE IF (data%State%MudSystem%MudType_MudElement%Array(i) == 4) THEN
            FlowEl(i)%MaterialType = 4
        ELSE
            FlowEl(i)%MaterialType = 1
        END IF
        
        IF (FlowEl(i)%MaterialType == 1 .OR. FlowEl(i)%MaterialType == 4) THEN
            FlowEl(i)%dPdLgrav = 0.0
        ELSE
            FlowEl(i)%dPdLgrav = 0.052 * FlowEl(i)%density
        END IF
        FlowEl(i)%StaticPressDiff = FlowEl(i)%dPdLgrav * REAL(FlowEl(i)%DepthDiff)


        
    END DO
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    !!!!!!!!!!!!!!! Well to choke manifold path
    j = 1
    DO i = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + data%State%FricPressDrop%NoAnnulusEl + 1 , data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + data%State%FricPressDrop%NoAnnulusEl + data%State%FricPressDrop%NoWelltoChokeEl
        FlowEl(i)%StartX = data%State%MudSystem%ChokeLine_Mud_Backhead_X%Array(j)
        FlowEl(i)%EndX = data%State%MudSystem%ChokeLine_Mud_Forehead_X%Array(j)
        FlowEl(i)%StartTVD = 0.d0
        FlowEl(i)%EndTVD = 0.d0
        FlowEl(i)%Id = 0.0
        FlowEl(i)%Od = data%Configuration%BopStack%ChokeLineId
        FlowEl(i)%density = data%State%MudSystem%ChokeLine_Density%Array(j)
            
        FlowEl(i)%Dhyd = FlowEl(i)%Od
        FlowEl(i)%Length = ABS(FlowEl(i)%EndX - FlowEl(i)%StartX)
        FlowEl(i)%DepthDiff = 0.d0
        FlowEl(i)%Area = PI / 4.0 * FlowEl(i)%Od**2 / Convfttoinch**2
        FlowEl(i)%Volume = FlowEl(i)%Area * REAL(FlowEl(i)%Length)
        IF (data%State%MudSystem%ChokeLine_MudOrKick%Array(j) == 0 .OR. data%State%MudSystem%ChokeLine_MudOrKick%Array(j) == 4) THEN      ! = 0 for mud, = 1 for gas kick, = 4 for air
            FlowEl(i)%MaterialType = data%State%MudSystem%ChokeLine_MudOrKick%Array(j)
        ELSE
            FlowEl(i)%MaterialType = 1
        END IF
        
            
        FlowEl(i)%dPdLgrav = 0.0
        FlowEl(i)%StaticPressDiff = 0.0
        
        !WRITE (*,*) 'Ch density, length, Type' , i, FlowEl(i)%density , FlowEl(i)%length, FlowEl(i)%MaterialType

        j = j + 1
            
    END DO
            
    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    !!!!!!!!!!!!!!! Open Hole
    data%State%FricPressDrop%OpenholeFirstEl = data%State%FricPressDrop%NoHorizontalEl + data%State%FricPressDrop%NoStringEl + data%State%FricPressDrop%NoAnnulusEl + data%State%FricPressDrop%NoWelltoChokeEl + 1
    j = 1
    DO i = data%State%FricPressDrop%OpenholeFirstEl , data%State%FricPressDrop%NumbEl ! = NoHorizontalEl + NoStringEl + NoCasingEl + NoOpenHoleEl
        FlowEl(i)%StartX = data%State%MudSystem%Xstart_OpMudElement%Array(j)
        FlowEl(i)%EndX = data%State%MudSystem%Xend_OpMudElement%Array(j)
        FlowEl(i)%StartTVD = data%State%MudSystem%TVDstart_OpMudElement%Array(j)
        FlowEl(i)%EndTVD = data%State%MudSystem%TVDend_OpMudElement%Array(j)
        FlowEl(i)%Id = 0.0
        FlowEl(i)%Od = data%State%MudSystem%PipeOD_OpMudElement%Array(j)
        FlowEl(i)%density = data%State%MudSystem%Density_OpMudElement%Array(j)
        FlowEl(i)%MaterialType = data%State%MudSystem%MudTypeOp_MudElement%Array(j)     ! = 0 for mud, = 1 for gas kick, = 4 for air
        
        FlowEl(i)%Dhyd = FlowEl(i)%Od
        FlowEl(i)%Length = ABS(FlowEl(i)%EndX - FlowEl(i)%StartX)
        FlowEl(i)%DepthDiff = ABS(FlowEl(i)%StartTVD - FlowEl(i)%EndTVD)
        FlowEl(i)%Area = PI / 4.0 * FlowEl(i)%Od**2 / Convfttoinch**2        
        FlowEl(i)%Volume = FlowEl(i)%Area * REAL(FlowEl(i)%Length)
        IF (data%State%MudSystem%MudTypeOp_MudElement%Array(j) == 0) THEN           ! = 0 for mud, = 1 for gas kick, = 4 for air
            FlowEl(i)%MaterialType = 0
        ELSE IF (data%State%MudSystem%MudTypeOp_MudElement%Array(j) == 4) THEN
            FlowEl(i)%MaterialType = 4
        ELSE
            FlowEl(i)%MaterialType = 1
        END IF
        
        j = j + 1
        IF (FlowEl(i)%MaterialType == 1) THEN
            FlowEl(i)%dPdLgrav = 0.0
        ELSE
            FlowEl(i)%dPdLgrav = 0.052 * FlowEl(i)%density
        END IF
        FlowEl(i)%StaticPressDiff = FlowEl(i)%dPdLgrav * REAL(FlowEl(i)%DepthDiff)
    
        !WRITE (*,*) 'O density, length, Type of ith element' , i, FlowEl(i)%density , FlowEl(i)%length, FlowEl(i)%MaterialType

    END DO
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
!!!!!!!!!!!!!!!!!!  Mud Element Viscosity Calculation
    
    DO i = 1 , data%State%FricPressDrop%NumbEl
        IF (FlowEl(i)%MaterialType /= 0)    CYCLE
        IF (FlowEl(i)%density < 8.33) THEN
            WRITE (*,*) 'H, S, A, Ch, O', data%State%FricPressDrop%NoHorizontalEl , data%State%FricPressDrop%StringFirstEl , data%State%FricPressDrop%AnnulusFirstEl , data%State%FricPressDrop%NoWellToChokeEl , data%State%FricPressDrop%OpenHoleFirstEl
            
            WRITE (*,*) 'Element No, Volume(gal), Density(ppg), Type'
            DO j = 1 , data%State%FricPressDrop%NumbEl
                WRITE (*,*) j, FlowEl(j)%Volume * ConvFt3toUSGal , FlowEl(j)%Density, FlowEl(j)%MaterialType
            END DO
        
            CALL ErrorStop('Density below water density (8.33)', i)
            CYCLE
        END IF
        
        FlowEl(i)%Theta600 = Theta600Refrence + 11.0 * (MAX(FlowEl(i)%density , 8.33) - DensityRefrence)
        FlowEl(i)%Theta300 = Theta300Refrence + 6.0 * (MAX(FlowEl(i)%density , 8.33) - DensityRefrence)
        !WRITE (*,*) 'i, Theta600, Theta600', i, FlowEl(i)%Theta600, FlowEl(i)%Theta300
    END DO
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    
    FlowEl(:)%FlowRate = 0.0
    FlowEl(:)%StartPress = 0.0
    FlowEl(:)%EndPress = 0.0
    
    !WRITE (*,*) 'NoGasPocket=' , SIZE(GasPocketWeight%Array)
    
    
    
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    data%State%MudSystem%OpenArea32 = ((data%EquipmentControl%ChokeManifold%LeftManualChoke * (1.0 - data%State%FricPressDrop%ManChoke1Plug * REAL(data%problems%ChokeProblems%ManualChoke1PluggedPercent)) / 100.0) + 0.5 * data%State%FricPressDrop%ManChoke1Washout) * data%State%Choke%ChokeAreaFullyOpen
    data%State%MudSystem%OpenArea33 = data%State%CHOOKE(1)%AreaChokeFinal * Convfttoinch**2              !  (1.d0 - data%State%CHOOKE(1)%PercentClose) * ChokeAreaFullyOpen
    data%State%MudSystem%OpenArea34 = data%State%CHOOKE(2)%AreaChokeFinal * Convfttoinch**2              !  (1.d0 - data%State%CHOOKE(2)%PercentClose) * ChokeAreaFullyOpen
    data%State%MudSystem%OpenArea35 = ((data%EquipmentControl%ChokeManifold%RightManualChoke * (1.0 - data%State%FricPressDrop%ManChoke2Plug * REAL(data%problems%ChokeProblems%ManualChoke2PluggedPercent)) / 100.0) + 0.5 * data%State%FricPressDrop%ManChoke2Washout) * data%State%Choke%ChokeAreaFullyOpen
    !write(*,*) 'OpenArea32=' , OpenArea32, active32, ManChoke1Plug, ManualChoke1PluggedPercent, ManChoke1Washout
    !write(*,*) 'OpenArea33=' , OpenArea33, active33
    !write(*,*) 'OpenArea34=' , OpenArea34, active34, !HydraulicChoke2PluggedPercent
    !write(*,*) 'OpenArea35=' , OpenArea35, active35, ManChoke2Plug, ManualChoke2PluggedPercent, ManChoke2Washout
    InstantaneousTotalOpenChokeArea = data%State%MudSystem%OpenArea32 * active32 + data%State%MudSystem%OpenArea33 * active33 + data%State%MudSystem%OpenArea34 * active34 + data%State%MudSystem%OpenArea35 * active35 + data%State%FricPressDrop%ChokeBypassArea * active29
    data%State%FricPressDrop%OldTotalOpenChokeArea = data%State%FricPressDrop%TotalOpenChokeArea
    !WRITE (*,*) 'Instantaneous / Old TotalOpenChokeArea', InstantaneousTotalOpenChokeArea, OldTotalOpenChokeArea
    IF (data%State%FricPressDrop%OldTotalOpenChokeArea <= 0.01 * data%State%Choke%ChokeAreaFullyOpen)     data%State%FricPressDrop%OldTotalOpenChokeArea = 0.01 * data%State%Choke%ChokeAreaFullyOpen
    IF (InstantaneousTotalOpenChokeArea <= 0.01 * data%State%Choke%ChokeAreaFullyOpen) THEN
        data%State%MudSystem%WellToChokeManifoldOpen = .FALSE.
        data%State%FricPressDrop%OldTotalOpenChokeArea = 0.01 * data%State%Choke%ChokeAreaFullyOpen
        !WRITE (*,*) ' Choke Controler Here 2'
        
    ELSE IF (InstantaneousTotalOpenChokeArea > 0.5 * data%State%Choke%ChokeAreaFullyOpen .OR. data%State%MudSystem%WellToPitsOpen) THEN
                    ! mud flows through well to bell nipple, or choke is rather open
        data%State%FricPressDrop%TotalOpenChokeArea = InstantaneousTotalOpenChokeArea
        !WRITE (*,*) ' Choke Controler Here 1'
    ELSE
        IF (InstantaneousTotalOpenChokeArea > data%State%FricPressDrop%OldTotalOpenChokeArea) THEN
            data%State%FricPressDrop%TotalOpenChokeArea = MIN(1.1 * data%State%FricPressDrop%OldTotalOpenChokeArea , InstantaneousTotalOpenChokeArea)
            WRITE (*,*) ' Choke is opening' !, TotalOpenChokeArea, OldTotalOpenChokeArea, InstantaneousTotalOpenChokeArea
        ELSE IF (InstantaneousTotalOpenChokeArea < data%State%FricPressDrop%OldTotalOpenChokeArea) THEN
            data%State%FricPressDrop%TotalOpenChokeArea = MAX(0.9 * data%State%FricPressDrop%OldTotalOpenChokeArea , InstantaneousTotalOpenChokeArea)
            WRITE (*,*) ' Choke is closing' !, TotalOpenChokeArea, OldTotalOpenChokeArea, InstantaneousTotalOpenChokeArea
        ELSE
            data%State%FricPressDrop%TotalOpenChokeArea = InstantaneousTotalOpenChokeArea
        END IF
    END IF
    !WRITE (*,*) 'TotalOpenChokeArea=', TotalOpenChokeArea 
    
        
    
    
            
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!            
    
!!!!!!!!!!!!!!!!!!!!  Well Head Condition ( Open or Closed )
    KickVARIABLES%ChokeKroneckerDelta = 0
    !WRITE (*,*) ' WelltoPitsOpen = ', WelltoPitsOpen
    !WRITE (*,*) ' WellToChokeManifoldOpen = ', WellToChokeManifoldOpen
    IF (data%State%MudSystem%WellToPitsOpen == .FALSE. .AND. data%State%MudSystem%WellToChokeManifoldOpen == .FALSE.) THEN
        !WRITE (*,*) ' WelltoPitsOpen = ', WelltoPitsOpen
        !WRITE (*,*) ' WellToChokeManifoldOpen = ', WellToChokeManifoldOpen
        
        KickVARIABLES%WellHeadOpen = .FALSE.
        !WRITE (*,*) ' Well Head is closed '
    ELSE
        KickVARIABLES%WellHeadOpen = .TRUE.
    END IF
    
    IF (data%State%MudSystem%WellToPitsOpen == .FALSE. .AND. data%State%MudSystem%WellToChokeManifoldOpen)        KickVARIABLES%ChokeKroneckerDelta = 1
    !WRITE (*,*) ' WelltoPitsOpen=' , WelltoPitsOpen, 'WellToChokeManifoldOpen= ' , WellToChokeManifoldOpen
    !WRITE (*,*) 'WellToChokeOpen , WellToChokeWasOpen', WellToChokeManifoldOpen, KickVARIABLES%WellToChokeManifoldWasOpen
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    IF (KickVARIABLES%NoGasPocket > 0) THEN
        IF (data%State%MudSystem%Choke_Kick_Saved_Volume_Final > 0.d0) THEN
            !WRITE (*,*) 'Choke_Kick_Saved_Volume_Final=', Choke_Kick_Saved_Volume_Final, (GasPocketModifiedVol%Array(NoGasPocket) * Convft3ToUSgal)
            KickVARIABLES%ExitMass = REAL(data%State%MudSystem%Choke_Kick_Saved_Volume_Final / (GasPocketModifiedVol%Array(KickVARIABLES%NoGasPocket) * Convft3ToUSgal)) * GasPocketWeight%Array(KickVARIABLES%NoGasPocket)
            !IF (NoGasPocket > 1)        ExitMass = ExitMass + (SUM(GasPocketDeltaVol%Array(1 : NoGasPocket - 1)) / GasPocketNewVol%Array(NoGasPocket)) * GasPocketWeight%Array(NoGasPocket) 
            !WRITE (*,*) 'Total Exit mass=' , ExitMass !, DeltaVolumePipe, data%State%MudSystem%StringFlowRate * dt / ConvMinToSec
            IF (FlowEl(data%State%FricPressDrop%OpenholeFirstEl - 1)%MaterialType /= 1) THEN
                CALL RemoveGasPocket(KickVARIABLES%NoGasPocket)
                WRITE (*,*) 'Choke_Kick_Saved_Volume_Final=', data%State%MudSystem%Choke_Kick_Saved_Volume_Final
            END IF
            
        END IF
        
        
        !NoGasPocket = SIZE(GasPocketWeight%Array)
        IF (KickVARIABLES%NoGasPocket > 0)         CALL GasPocketFlowElementTransformer
    END IF
    
    IF (data%State%FricPressDrop%NoWelltoChokeEl > 0) THEN
        KickVARIABLES%OldChokeDensity = KickVARIABLES%ChokeDensity
        VolumeDensityProduct = 0.0
        Do i = data%State%FricPressDrop%AnnulusLastEl + 1 , data%State%FricPressDrop%OpenholeFirstEl - 1
            VolumeDensityProduct = VolumeDensityProduct + ((MAX(FlowEl(i)%Density , KickVARIABLES%ChokeMinDensity)) * FlowEl(i)%Volume)
        End Do
        KickVARIABLES%ChokeDensity = VolumeDensityProduct / SUM(FlowEl(data%State%FricPressDrop%AnnulusLastEl + 1 : data%State%FricPressDrop%OpenholeFirstEl - 1)%Volume)
                            
        !WRITE (*,*) 'ChokeDensity=', KickVARIABLES%ChokeDensity
        data%State%FricPressDrop%Kchoke = (KickVARIABLES%ChokeDensity / ((2.0 * 89158.0) * (0.26 * 0.61 * data%State%FricPressDrop%TotalOpenChokeArea)**2)) * 4.d0    ! *4.d0: seyyed gofte
    END IF
    
    IF (data%State%FricPressDrop%BitTotallyPluged == .FALSE.) THEN
        data%State%FricPressDrop%KBit = FlowEl(data%State%FricPressDrop%StringLastEl)%density / 12042.0 / BitCd**2 / BitTotNozzArea**2
    END IF
    
    
    
END SUBROUTINE WellPressureDataTransfer
    
    
SUBROUTINE GasPocketFlowElementTransformer

    !! This Subroutine makes relationship between fluid flow elements and gas pocket elements
    !! The main variable of this subroutine is GasPocketFlowEl
    !! For each gas pocket, GasPocketFlowEl has a row, which columns are numbers of flow elements which are contains that gas pocket
    !! For example, if we have two gas pockets and they are in (7,8) and (12,3,4) flow elements respectively,
    !! GasPocketFlowEl(1,:) = [12 , 3 , 4] and GasPocketFlowEl(2,:) = [7 , 8 , 0]


    USE FricPressDropVarsModule
    USE MudSystemVARIABLES
use SimulationVariables !@@@
    USE GeoElements_FluidModule
    USE Fluid_Flow_Startup_Vars
    use KickVARIABLESModule
    USE CMudPropertiesVariables
    USE CBopStackVariables
    USE CDownHoleVariables
    USE CError
    
    IMPLICIT NONE
    
    INTEGER :: i , j , k
    REAL :: PressureCorrection
    PressureCorrection = 1.0
    
    ALLOCATE(KickVARIABLES%GasPocketFlowEl(KickVARIABLES%NoGasPocket , 1))
    KickVARIABLES%GasPocketFlowEl(:,:) = 0
    !WRITE (*,*) 'NoGasPocket=' , NoGasPocket
    
    i = 1
    j = 1
    DO k = data%State%FricPressDrop%OpenholeFirstEl , data%State%FricPressDrop%NumbEl
        IF (FlowEl(k)%MaterialType == 1) THEN
            IF (j > 1) THEN
                IF (KickVARIABLES%GasPocketFlowEl(i , j - 1) /= k - 1) THEN
                    i = i + 1
                    j = 1
                END IF
            END IF
            
            IF (SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2) < j) THEN
                ALLOCATE(KickVARIABLES%tempGasPocketFlowEl(SIZE(KickVARIABLES%GasPocketFlowEl , dim = 1) , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)))
                KickVARIABLES%tempGasPocketFlowEl = KickVARIABLES%GasPocketFlowEl
                DEALLOCATE(KickVARIABLES%GasPocketFlowEl)
                ALLOCATE(KickVARIABLES%GasPocketFlowEl(SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 1) , SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 2) + 1))
                KickVARIABLES%GasPocketFlowEl(: , 1 : SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 2)) = KickVARIABLES%tempGasPocketFlowEl
                KickVARIABLES%GasPocketFlowEl(: , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)) = 0
                DEALLOCATE(KickVARIABLES%tempGasPocketFlowEl)
            END IF
            
            KickVARIABLES%GasPocketFlowEl(i , j) = k
            j = j + 1
        END IF
    END DO
    !WRITE (*,*) 'GasPocketFlowEl=' , GasPocketFlowEl
    
    DO k = data%State%FricPressDrop%AnnulusFirstEl , data%State%FricPressDrop%OpenholeFirstEl - 1
        IF (FlowEl(k)%MaterialType == 1) THEN
            IF (j > 1) THEN
                IF (k == data%State%FricPressDrop%AnnulusFirstEl  .AND.  KickVARIABLES%GasPocketFlowEl(i , j - 1) == data%State%FricPressDrop%NumbEl) THEN
                    !WRITE (*,*) 'Kick is Around Bit'
                ELSE IF (KickVARIABLES%GasPocketFlowEl(i , j - 1) /= k - 1) THEN
                    i = i + 1
                    j = 1
                    !WRITE (*,*) 'i, j, k', i, j, k, FlowEl(k)%MaterialType
                END IF
            END IF
            
            IF (j > SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)) THEN
                ALLOCATE(KickVARIABLES%tempGasPocketFlowEl(SIZE(KickVARIABLES%GasPocketFlowEl , dim = 1) , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)))
                KickVARIABLES%tempGasPocketFlowEl = KickVARIABLES%GasPocketFlowEl
                DEALLOCATE(KickVARIABLES%GasPocketFlowEl)
                ALLOCATE(KickVARIABLES%GasPocketFlowEl(SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 1) , SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 2) + 1))
                KickVARIABLES%GasPocketFlowEl(: , 1 : SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 2)) = KickVARIABLES%tempGasPocketFlowEl
                KickVARIABLES%GasPocketFlowEl(: , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)) = 0
                DEALLOCATE(KickVARIABLES%tempGasPocketFlowEl)
            END IF
            IF (i > KickVARIABLES%NoGasPocket) THEN
                WRITE (*,*) 'GasPocketFlowEl', KickVARIABLES%GasPocketFlowEl
                WRITE (*,*) 'i , j , k', i, j, k
                WRITE (*,*) 'H, S, A, Ch, O', data%State%FricPressDrop%NoHorizontalEl , data%State%FricPressDrop%StringFirstEl , data%State%FricPressDrop%AnnulusFirstEl , data%State%FricPressDrop%NoWellToChokeEl , data%State%FricPressDrop%OpenHoleFirstEl
            
                WRITE (*,*) 'Element No, Volume(gal), Density(ppg),         Type'
                DO j = 1 , data%State%FricPressDrop%NumbEl
                    WRITE (*,*) j, FlowEl(j)%Volume * ConvFt3toUSGal , FlowEl(j)%Density, FlowEl(j)%MaterialType
                END DO

                CALL ErrorStop('Error in calculating GasPocketFlowEl')
            END IF

            KickVARIABLES%GasPocketFlowEl(i , j) = k
            !WRITE (*,*) 'GasPocketFlowEl=' , GasPocketFlowEl
            j = j + 1
            
        END IF
        !IF (i > NoGasPocket)        CALL ErrorStop('Error in GasPocketFlowEl' , i)
    END DO
    

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!    
    
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!    
    GasPocketModifiedVol%Array(:) = 0.d0
    
    DO i = 1 , KickVARIABLES%NoGasPocket
        !WRITE(*,*) 'GasPocketFlowEl Data' , i , GasPocketFlowEl(: , i)
        DO j = 1 , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)
            IF (KickVARIABLES%GasPocketFlowEl(i , j) > 0)     GasPocketModifiedVol%Array(i) = GasPocketModifiedVol%Array(i) + FlowEl(KickVARIABLES%GasPocketFlowEl(i , j))%Volume
        END DO
    END DO
    i = data%State%FricPressDrop%AnnulusLastEl
    j = data%State%FricPressDrop%OpenholeFirstEl - 1
    
    
    GasPocketOldVol%Array(:) = GasPocketModifiedVol%Array(:)
    
    GasPocketOldPress%Array(:) = GasPocketNewPress%Array(:)
    !WRITE (*,*) ' GasPocketweight%Array(1)', GasPocketweight%Array(1)
    IF (GasPocketWeight%Array(1) < 0.0)     CALL ErrorStop(' Error in Calculating Mass ' , KickVARIABLES%KickmdotBCoef)
    
    
    GasPocketDensity%Array(:) = (GasPocketweight%Array(:) / GasPocketModifiedVol%Array(:)) / convft3toUSgal     ! [lbm/ft^3 to ppg]
    DO i = 1 , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 1)
        DO j = 1 , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)
            IF (KickVARIABLES%GasPocketFlowEl(i , j) > 0)     FlowEl(KickVARIABLES%GasPocketFlowEl(i , j))%Density = GasPocketDensity%Array(i)
        END DO
    END DO
    
    
    !WRITE (*,*) 'Kick density (ppg)=' , GasPocketDensity(1)
    DownHole%InfluxRate = MAX(((KickVARIABLES%KickmdotACoef * (KickVARIABLES%KickmdotBCoef - GasPocketNewPress%Array(1))) / GasPocketDensity%Array(1) * ConvMinToSec) , 0.0)
    !WRITE (*,*) ' InfluxRate (gpm) =', InfluxRate
    
    i = data%State%FricPressDrop%OpenholeFirstEl - 1
    IF (data%State%MudSystem%Choke_Kick_Saved_Volume_Final > 0.d0 .AND. KickVARIABLES%WellHeadOpen) THEN
        WRITE (*,*) 'Kick is exiting through the choke'
        !KickVARIABLES%ExitMass = GasPocketDensity%Array(NoGasPocket) * (SUM(GasPocketFlowInduced%Array(:)) + (DeltaVolumePipe * ConvMinToSec / dt) + data%State%MudSystem%StringFlowRate) / ConvMinToSec * dt
        GasPocketWeight%Array(KickVARIABLES%NoGasPocket) = GasPocketweight%Array(KickVARIABLES%NoGasPocket) - KickVARIABLES%ExitMass
        IF (GasPocketWeight%Array(KickVARIABLES%NoGasPocket) > 0.0) THEN
            GasPocketOldPress%Array(KickVARIABLES%NoGasPocket) = GasPocketWeight%Array(KickVARIABLES%NoGasPocket) * data%State%GasType(KickVARIABLES%KickGasType)%GasConstant * &
                                                    GasPocketCompressibility%Array(KickVARIABLES%NoGasPocket) * GasPocketNewTemp%Array(KickVARIABLES%NoGasPocket) / GasPocketModifiedVol%Array(KickVARIABLES%NoGasPocket)
            GasPocketDensity%Array(KickVARIABLES%NoGasPocket) = (GasPocketweight%Array(KickVARIABLES%NoGasPocket) / GasPocketModifiedVol%Array(KickVARIABLES%NoGasPocket)) / convft3toUSgal
            GasPocketNewVol%Array(:) = GasPocketOldVol%Array(:)
        ELSE    ! gas pocket is escaped from the well completely
            WRITE (*,*) ' Last Pocket Removed'
            WRITE (*,*) ' GasPocketFlowEl', KickVARIABLES%GasPocketFlowEl
            CALL RemoveGasPocket(KickVARIABLES%NoGasPocket)
            
            !WRITE (*,*) ' GasPocketFlowEl', GasPocketFlowEl
            
        END IF
        
    ELSE IF (KickVARIABLES%NoGasPocket == 1 .OR. KickVARIABLES%WellHeadOpen) THEN        ! kick is not last element of choke line or well head is closed
            GasPocketNewPress%Array(:) = GasPocketNewPress%Array(:) * REAL(GasPocketNewVol%Array(:) / GasPocketOldVol%Array(:))
            GasPocketOldPress%Array(:) = GasPocketNewPress%Array(:)
            GasPocketNewVol%Array(:) = GasPocketOldVol%Array(:)
    
    ELSE    ! Many gas pocket is in the well and well head is closed
        
103     FORMAT (2I, 4X, (F9.3), 5X, (F9.3))
        
        GasPocketNewPress%Array(KickVARIABLES%NoGasPocket) = GasPocketOldPress%Array(KickVARIABLES%NoGasPocket)
        GasPocketNewVol%Array(KickVARIABLES%NoGasPocket) = GasPocketCompressibility%Array(KickVARIABLES%NoGasPocket) * data%State%GasType(KickVARIABLES%KickGasType)%GasConstant * &
                            GasPocketNewTemp%Array(KickVARIABLES%NoGasPocket) * GasPocketWeight%Array(KickVARIABLES%NoGasPocket) / GasPocketNewPress%Array(KickVARIABLES%NoGasPocket)
        
        
        
                !WRITE (*,*) 'Mid'
        DO WHILE (ABS(PressureCorrection * 10.0) > KickConvergenceTolerance)
            DO j = KickVARIABLES%NoGasPocket - 1 , 1 , -1
                CALL KickFunctionsCalculator(GasPocketNewPress%Array(j) , j , 5)
                GasPocketNewVol%Array(j) = GasPocketCompressibility%Array(j) * data%State%GasType(KickVARIABLES%KickGasType)%GasConstant * &
                            GasPocketNewTemp%Array(j) * GasPocketWeight%Array(j) / GasPocketNewPress%Array(j)
                !WRITE (*,*) j, REAL(GasPocketNewPress%Array(j)), REAL(GasPocketNewVol%Array(j) * Convft3toUSgal)
            END DO
            PressureCorrection = (SUM(GasPocketNewVol%Array(:)) - SUM(GasPocketOldVol%Array(:))) / SUM(GasPocketNewVol%Array(:) / GasPocketNewPress%Array(:))
            GasPocketNewPress%Array(KickVARIABLES%NoGasPocket) = GasPocketNewPress%Array(KickVARIABLES%NoGasPocket) + PressureCorrection
            GasPocketNewVol%Array(KickVARIABLES%NoGasPocket) = GasPocketCompressibility%Array(KickVARIABLES%NoGasPocket) * data%State%GasType(KickVARIABLES%KickGasType)%GasConstant * &
                            GasPocketNewTemp%Array(KickVARIABLES%NoGasPocket) * GasPocketWeight%Array(KickVARIABLES%NoGasPocket) / GasPocketNewPress%Array(KickVARIABLES%NoGasPocket)
            j = KickVARIABLES%NoGasPocket
            !WRITE (*,*) j, REAL(GasPocketNewPress%Array(j)), REAL(GasPocketNewVol%Array(j) * Convft3toUSgal)
        END DO
        
        !WRITE (*,*) 'Pressure Correction = ', PressureCorrection
        WRITE (*,*) 'After'
        DO i = 1 , KickVARIABLES%NoGasPocket
            WRITE (*,*) i, REAL(GasPocketNewPress%Array(i)), REAL((GasPocketNewVol%Array(i) * Convft3toUSgal))
        END DO
        
        GasPocketOldPress%Array(:) = GasPocketNewPress%Array(:)
        
    END IF
    
    
    IF (KickVARIABLES%GasPocketFlowEl(KickVARIABLES%NoGasPocket , 1) == 0) THEN
        CALL RemoveGasPocket(1)
        
        !WRITE (*,*) 'GasPocketFlowEl', GasPocketFlowEl
        WRITE (*,*) 'First Pocket Removed'
        !WRITE (*,*) 'GasPocketFlowEl', GasPocketFlowEl
    END IF
    
    
    
    DownHole%KickVolume = INT(SUM(GasPocketOldVol%Array(:)) * convft3toUSgal / 42. * 10.0) / 10.0
    
    
    !WRITE (*,*) ' Gas Kick Pressure (psi) = ' , GasPocketOldPress(1) , INT((GasPocketNewVol(1) / GasPocketOldVol(1)) * 1000.d0) / 1000.d0
    !WRITE (*,*) ' b) GasPocketNewVol (ft^3) = ' , GasPocketOldVol(1)
        
        
    
    END SUBROUTINE
    
    
SUBROUTINE RemoveGasPocket(ilocal)
    
    use KickVARIABLESModule
    USE MudSystemVARIABLES
use SimulationVariables !@@@
    
    IMPLICIT NONE
    
    INTEGER :: ilocal

    DEALLOCATE(KickVARIABLES%KickJacobian , KickVARIABLES%OldKickJacobian , KickVARIABLES%KickVandPFunction , KickVARIABLES%KickUnknownVector , KickVARIABLES%KickCorrectionVector)
    
        
    WRITE (*,*) 'Gas Pocket Will be Removed', ilocal
    call GasPocketOldPress%Remove(ilocal)
    call GasPocketNewPress%Remove(ilocal)
    call GasPocketOldTemp%Remove(ilocal)
    call GasPocketNewTemp%Remove(ilocal)
    call GasPocketOldVol%Remove(ilocal)
    call GasPocketNewVol%Remove(ilocal)
    call GasPocketDeltaVol%Remove(ilocal)
    call GasPocketFlowInduced%Remove(ilocal)
    call GasPocketModifiedVol%Remove(ilocal)
    call GasPocketWeight%Remove(ilocal)
    call GasPocketDensity%Remove(ilocal)
    call GasPocketCompressibility%Remove(ilocal)
    
    KickVARIABLES%NoGasPocket = KickVARIABLES%NoGasPocket - 1
     
    IF (KickVARIABLES%NoGasPocket > 0) THEN
        ALLOCATE(KickVARIABLES%KickJacobian(2 * KickVARIABLES%NoGasPocket , 2 * KickVARIABLES%NoGasPocket) , KickVARIABLES%OldKickJacobian(2 * KickVARIABLES%NoGasPocket , 2 * KickVARIABLES%NoGasPocket) , KickVARIABLES%KickVandPFunction(2 * KickVARIABLES%NoGasPocket))
        ALLOCATE(KickVARIABLES%KickUnknownVector(2 * KickVARIABLES%NoGasPocket) , KickVARIABLES%KickCorrectionVector(2 * KickVARIABLES%NoGasPocket))
        
        IF (ALLOCATED(KickVARIABLES%GasPocketFlowEl)) THEN
            ALLOCATE(KickVARIABLES%tempGasPocketFlowEl(KickVARIABLES%NoGasPocket , SIZE(KickVARIABLES%GasPocketFlowEl , dim = 2)))
            KickVARIABLES%tempGasPocketFlowEl = KickVARIABLES%GasPocketFlowEl(1 : KickVARIABLES%NoGasPocket , :)
            DEALLOCATE(KickVARIABLES%GasPocketFlowEl)
            ALLOCATE(KickVARIABLES%GasPocketFlowEl(SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 1) , SIZE(KickVARIABLES%tempGasPocketFlowEl , dim = 2)))
            KickVARIABLES%GasPocketFlowEl = KickVARIABLES%tempGasPocketFlowEl
            DEALLOCATE(KickVARIABLES%tempGasPocketFlowEl)
        END IF
        
    ELSE        ! NoGasPocket = 0
        
        
        CALL GasPocketOldPress%Empty
        CALL GasPocketNewPress%Empty
        CALL GasPocketOldTemp%Empty
        CALL GasPocketNewTemp%Empty
        CALL GasPocketOldVol%Empty
        CALL GasPocketNewVol%Empty
        CALL GasPocketDeltaVol%Empty
        CALL GasPocketFlowInduced%Empty
        CALL GasPocketModifiedVol%Empty
        CALL GasPocketWeight%Empty
        CALL GasPocketDensity%Empty
        
        DownHole%InfluxRate = 0.0
        DownHole%KickVolume = 0.0
                
    END IF         
    
    !WRITE (*,*) 'Gas Pocket Removed', NoGasPocket + 1, ALLOCATED(GasPocketWeight%Array)
    
END SUBROUTINE