SimulationCore2/Rop/ROP_MainCalculation.f90

subroutine ROP_MainCalculation
    use ConfigurationVariables
    use sROP_Other_Variables
    use sROP_Variables
    use CStringConfigurationVariables
    use CformationVariables
    ! use CSimulationVariables
    use CDataDisplayConsoleVariables
    use CDrillingConsoleVariables
    use CmudPropertiesVariables
    use CHoistingVariables
    ! use CDrillingConsole
    use CPathGenerationVariables
    use RTable_VARIABLES, only: RTable
    use TD_DrillStemComponents
    use TD_WellGeometry
    use PressureDisplayVARIABLESModule
    use MudSystemVARIABLES
    USE FricPressDropVarsModule
    use CReservoirVariables
    use CWarningsVariables
    use TopDrive_VARIABLES, only: TDS
    use TD_GeneralData
    
    implicit none
    
    Integer :: i , zero_ROPcount
    !Real(8) :: ROP_Bit%SetROPGauge
    
    zero_ROPcount = 0
    ROP_Bit%NoOfFormations = Configuration%Formation%Count
    ROP_Bit%DrillingVerticalDepth = TD_WellGeneral%WellTotalVerticalLength
    
    
    
   !===> MaximumWellDepthExceeded Warning
    if ( ROP_Bit%DrillingVerticalDepth>=(Configuration%Formation%Formations(Configuration%Formation%Count)%Top+Configuration%Formation%Formations(Configuration%Formation%Count)%Thickness) ) then
        ROP_Bit%RateOfPenetration = 0.0d0
        Call Set_ROP(ROP_Bit%RateOfPenetration)
        Call Activate_MaximumWellDepthExceeded()
        return
    end if
   !=====================================
    
    
    
    if ( ROP_Spec%FormationNumber/=0 .and. DrillingConsole%HideDrillingBrake==1 ) then      ! Hide Drilling Brake Mode
        ROP_Spec%FormationNumber = ROP_Spec%FormationNumber
    else
        do i= 1,ROP_Bit%NoOfFormations   
            ROP_Spec%FormationTopDepth = Configuration%Formation%Formations(i)%Top
            if (ROP_Bit%DrillingVerticalDepth>=ROP_Spec%FormationTopDepth) then    
                ROP_Spec%FormationNumber = i   
            end if
        end do
    end if
    
    
    
    
    !!===> Hide Drilling Brake Mode
    !if ( ROP_Spec%FormationNumber==FormationNo .and. HideDrillingBrake==1 ) then  !????????????
    !    ROP_Spec%FormationNumber = FormationNo-1
    !end if
    !!=============================
    
    
    
    
    !ROP_Spec%BitClass = BitDefinition%BitCode     !????????????
    call Bit_Specification
    
    
    
    
    ! $$$$$**$$$$$**$$$$$**$$$$$**$$$$$**     Variables Initialization:    *$$$$$**$$$$$**$$$$$**$$$$$**$$$$$  
    ROP_Spec%DiameterOfBit        =  Configuration%StringConfiguration%BitDefinition%BitSize          ! unit : [in.] (Typical Range: 3.0 to 30.0)
    ROP_Spec%NumberOfBitNozzles   =  Configuration%StringConfiguration%BitDefinition%BitNozzleNo      ! (Typical Values: 1 to 10)
    ROP_Spec%DiameterOfBitNozzle  =  Configuration%StringConfiguration%BitDefinition%BitNozzleSize    ! unit : [inch] *** basic input: [1/32 in.] (Typical Range: 8.0 to 32.0)
    ROP_Spec%CriticalMudDensity   =  Configuration%Formation%Formations(ROP_Spec%FormationNumber)%PorePressureGradient/.465d0*9.d0      ! ????????? delete ,unit : [ppg] or [lb/gal] (Typical Range: 0 to 10.0)
    ROP_Bit%FormationMudDensity   =  Configuration%Formation%Formations(ROP_Spec%FormationNumber)%PorePressureGradient/0.052d0
    ROP_Bit%BottomHolePressure    =  PressureDisplayVARIABLES%PressureGauges(3)              !5200 [psi]
    ROP_Bit%ECD                   =  ROP_Bit%BottomHolePressure/(0.052*ROP_Bit%DrillingVerticalDepth)
    ROP_Spec%CriticalWeightOnBit  =  (Configuration%Formation%Formations(ROP_Spec%FormationNumber)%ThresholdWeight/5.d0)-(.06d0*(Configuration%Formation%Formations(ROP_Spec%FormationNumber)%ThresholdWeight-10.d0))      ! unit : [klb/in] (Typical Range: 0 to 10 ----> 0.6 to 2)
    !IF (ALLOCATED(FlowEl)) THEN
    !    ROP_Bit%MudViscosity     =  FlowEl(NoHorizontalEl + NoStringEl)%mueff     !13.5  [cP]
        ROP_Bit%MudDensity        = MudSystem%BitMudDensity ! [ppg]
    !ELSE
        ROP_Bit%MudViscosity      =  13.5         ! [cP]
        !ROP_Bit%MudDensity       =  9.2          ! [ppg]
    !END IF
    ROP_Bit%MudFlowrate           =  MudSystem%StringFlowRateFinal    ! [gpm]
    ROP_Spec%ReynoldsNumber        =  ROP_Bit%MudFlowrate*ROP_Bit%MudDensity/(ROP_Bit%MudViscosity*ROP_Spec%NumberOfBitNozzles*ROP_Spec%DiameterOfBitNozzle)           ! unit : [dimensionless] (Typical Range: 0.1 to 1000.0)    
    ! $$$$$**$$$$$**$$$$$**$$$$$**$$$$    End of Variable Initialization    $$$$$**$$$$$**$$$$$**$$$$$**$$$$$
    
    
    
    
    
    ! -----**-----**-----**-----**-----*      ROP_bit%RateofPenetration Model Coefficients:     *-----**-----**-----**-----**-----
    ROP_Spec%a1 = log(Configuration%Formation%Formations(ROP_Spec%FormationNumber)%Drillablity) 
    ROP_Spec%a2 = 1.2799d-04
    ROP_Spec%a3 = 1.7952d-04
    ROP_Spec%a4 = 4.0656d-05
    ROP_Spec%a5 = 2.9021d-01
    ROP_Spec%a6 = 9.4882d-02
    ROP_Spec%a7 = 2.1837d-01
    ROP_Spec%a8 = 4.4915d-01
    ROP_Spec%dt = 0.1d0   ![s]
    ROP_Spec%TouH = Configuration%Formation%Formations(ROP_Spec%FormationNumber)%Abrasiveness*3600.d0    ! [hr]--->[s]   ( Typical Range: 1<Abrasiveness<100 )
    ! -----**-----**-----**-----**---     End of ROP_bit%RateofPenetration Model Coefficients      ---**-----**-----**-----**-----
    
    
    
    
    
    
    ! $$$$$**$$$$$**$$$$$**$$$$$**$$$$$**        Main Calculations:        *$$$$$**$$$$$**$$$$$**$$$$$**$$$$$
    ROP_Spec%f1 = Configuration%Formation%Formations(ROP_Spec%FormationNumber)%Drillablity                        ! 1<=Drillability<=100  [ft/h]
    ROP_Spec%f2 = exp(2.303d0*ROP_Spec%a2*(10000.d0 - ROP_Bit%DrillingVerticalDepth))            ! First Compaction Vairable
    ROP_Spec%f3 = exp(2.303d0*ROP_Spec%a3*(ROP_Bit%DrillingVerticalDepth**0.69d0)*(ROP_Bit%FormationMudDensity-9.d0))    ! Second Compaction Vairable
    
    
    if ( (2.303d0*ROP_Spec%a4*ROP_Bit%DrillingVerticalDepth*(ROP_Bit%FormationMudDensity-ROP_Bit%ECD))>=0.d0 ) then
        ROP_Spec%f4 = exp(2.303d0*ROP_Spec%a4*ROP_Bit%DrillingVerticalDepth*(ROP_Bit%FormationMudDensity-ROP_Bit%ECD))           ! Underbalance Drilling Variable
    else
        ROP_Spec%f4 = 1.0d0
    end if
    
    
    if (TD_String%WeightOnBit>0.d0) then
        ROP_Bit%WeightOnBit = TD_String%WeightOnBit/1000.d0    ![klb]
    else
        ROP_Bit%WeightOnBit = 0.d0
    end if
    
    if ( (ROP_Bit%WeightOnBit/ROP_Spec%DiameterOfBit)<ROP_Spec%CriticalWeightOnBit ) then
        ROP_Spec%f5 = 0.0d0
    else
        ROP_Spec%f5 = ( ((ROP_Bit%WeightOnBit/ROP_Spec%DiameterOfBit)-ROP_Spec%CriticalWeightOnBit)/(4.d0-ROP_Spec%CriticalWeightOnBit) )**ROP_Spec%a5 
    end if
   
    
    if (Hoisting%DriveType==1 .and. RTable%Speed>0.d0) then
        ROP_Bit%RotarySpeed = RTable%Speed               ![rpm]
    else if (Hoisting%DriveType==0 .and. (TDS%Speed>0. .or. RTable%Speed>0.)) then
        ROP_Bit%RotarySpeed = TDS%Speed+RTable%Speed     ![rpm]
    else
        ROP_Bit%RotarySpeed = 0.0d0
    end if
    ROP_Spec%f6 = (ROP_Bit%RotarySpeed/100.d0)**ROP_Spec%a6
    
    
    ROP_Spec%f7 = exp(-ROP_Spec%a7*ROP_Bit%BitWearing)
    
    
    Call JetImpactForce
    ROP_Spec%f8 = ROP_Bit%JetImpactForce/1000.d0
    
    
    
    
    ROP_Bit%RateOfPenetration = (ROP_Spec%f1*ROP_Spec%f2*ROP_Spec%f3*ROP_Spec%f4*ROP_Spec%f5*ROP_Spec%f6*ROP_Spec%f7*ROP_Spec%f8)    ![ft/h]
    ROP_Bit%RateOfPenetration = (DINT(ROP_Bit%RateOfPenetration*10.d0))/10.d0
    
    if ( (TD_WellGeneral%WellTotalLength==PathGeneration%Items(PathGeneration%ItemCount)%MeasuredDepth) ) then
        ROP_Bit%SetROPGauge = ROP_Bit%RateOfPenetration
        Call Set_ROP(ROP_Bit%SetROPGauge)                     ![ft/h]
        ROP_Bit%OldROPValue(4) = ROP_Bit%RateOfPenetration
        !print* , 'first rop=' , ROP_Bit%OldROPValue , ROP_Bit%RateOfPenetration ,&
        !    zero_ROPcount , ROP_Bit%SetROPGauge , ROP_Bit%OldROPDepth , TD_WellGeneral%WellTotalLength , PathGenerations(PathGenerationCount)%MeasuredDepth
    else if ( ((TD_WellGeneral%WellTotalLength+(ROP_Bit%RateOfPenetration*TD_General%TimeStep/3600.d0))-ROP_Bit%OldROPDepth)>=0.1 ) then
        do i= 1,4
            if ( ROP_Bit%OldROPValue(i)==0. ) then
                zero_ROPcount = zero_ROPcount+1
            end if
        end do
        ROP_Bit%SetROPGauge = (ROP_Bit%RateOfPenetration+ROP_Bit%OldROPValue(1)+ROP_Bit%OldROPValue(2)+ROP_Bit%OldROPValue(3)+ROP_Bit%OldROPValue(4))/sngl(5-zero_ROPcount)
        Call Set_ROP(ROP_Bit%SetROPGauge)                     ![ft/h]
        do i= 2,4
            ROP_Bit%OldROPValue(i-1) = ROP_Bit%OldROPValue(i)
        end do
        ROP_Bit%OldROPValue(4) = ROP_Bit%RateOfPenetration
        ROP_Bit%OldROPDepth    = TD_WellGeneral%WellTotalLength+(ROP_Bit%RateOfPenetration*TD_General%TimeStep/3600.d0)
        !print* , 'new rop=' , ROP_Bit%OldROPValue , ROP_Bit%RateOfPenetration ,&
        !    zero_ROPcount , ROP_Bit%SetROPGauge , ROP_Bit%OldROPDepth , TD_WellGeneral%WellTotalLength , PathGenerations(PathGenerationCount)%MeasuredDepth
    else
        Call Set_ROP(ROP_Bit%SetROPGauge)                     ![ft/h]
        !print* , 'old rop=' , ROP_Bit%OldROPValue , ROP_Bit%RateOfPenetration ,&
        !    zero_ROPcount , ROP_Bit%SetROPGauge , ROP_Bit%OldROPDepth , TD_WellGeneral%WellTotalLength , PathGenerations(PathGenerationCount)%MeasuredDepth
    end if    
    
    
    
    
    
    
    if (ROP_Bit%RotarySpeed > 0.d0) THEN 
        ROP_Bit%BitTorque = ( 3.79d0 + 19.17d0*sqrt( ROP_Bit%RateOfPenetration / (ROP_Bit%RotarySpeed*ROP_Spec%DiameterOfBit)) ) * ROP_Spec%DiameterOfBit * ROP_Bit%WeightOnBit * ( 1.d0 / ( 1.d0 + 0.00021d0*ROP_Bit%DrillingVerticalDepth) )
        !ROP_Bit%BitTorque = ROP_Bit%BitTorque/3.   !bi dalil taghsim bar 3 shode(chon adad bozorg bude), baadan az rabete check shavad (seyyed gofte)
    else
        ROP_Bit%BitTorque = 0.d0
    end if
   
    
    
    if ( (ROP_Bit%WeightOnBit/ROP_Spec%DiameterOfBit)<(ROP_Spec%wdmax) ) then
        ROP_Bit%BitWearing = ROP_Bit%BitWearing +( (ROP_Spec%dt*ROP_Spec%H3/ROP_Spec%TouH)*((ROP_Bit%RotarySpeed/100.d0)**ROP_Spec%H1)*((ROP_Spec%wdmax-4.d0)/(ROP_Spec%wdmax-(ROP_Bit%WeightOnBit/ROP_Spec%DiameterOfBit)))*((1.d0+(ROP_Spec%H2/2.d0))/(1.d0+(ROP_Spec%H2*ROP_Bit%BitWearing))) )
    else
        ROP_Bit%BitWearing = 1.0d0          !( Typical Range: 0<=ROP_Bit%BitWearing<=1 )
    end if
    
    
    
    ROP_Bit%BearingWear = ROP_Bit%BearingWear+(ROP_Spec%dt/3600.d0)*(ROP_Bit%RotarySpeed/100.d0/ROP_Spec%BrCoef)*((ROP_Bit%WeightOnBit/4.d0/ROP_Spec%DiameterOfBit)**1.5d0)
    
    
    !print*, 'ROP_Bit%RateOfPenetration=', ROP_Bit%RateOfPenetration
    !!print*, 'ROP_Bit%FormationMudDensity=', ROP_Bit%FormationMudDensity
    !!print*, 'ROP_Bit%ECD=', ROP_Bit%ECD
    !!print*, 'ROP_Bit%DrillingVerticalDepth=', ROP_Bit%DrillingVerticalDepth
    !!print*, 'power=', (2.303*ROP_Spec%a4*ROP_Bit%DrillingVerticalDepth*(ROP_Bit%FormationMudDensity-ROP_Bit%ECD))
    !print*, 'ROP_Bit%RotarySpeed=', ROP_Bit%RotarySpeed
    !!
    !print*, 'ROP_Spec%f1=', ROP_Spec%f1
    !print*, 'ROP_Spec%f2=', ROP_Spec%f2
    !print*, 'ROP_Spec%f3=', ROP_Spec%f3
    !print*, 'ROP_Spec%f4=', ROP_Spec%f4
    !print*, 'ROP_Spec%f5=', ROP_Spec%f5
    !print*, 'ROP_Spec%f6=', ROP_Spec%f6
    !print*, 'ROP_Spec%f7=', ROP_Spec%f7
    !print*, 'ROP_Spec%f8=', ROP_Spec%f8
    !print*, '***********************'

    
    
    
    
end subroutine