New Math Functions in Fortran

Posted
Comments None

For the Fortran code of the 39 new math functions added to the built-in math functions of GeneXproTools 5.0 with Mini-Release 1, a good starting point is the Visual Basic code. Although there are quite a few differences between Fortran and Visual Basic, the basic structures that we need for the Fortran code are already there, namely, IF THEN ELSE structures with END IF and using the function name to return a function output.

Another similarity between both languages concerns the initialization of variables: in both cases variables cannot be initialized during declaration. But in Fortran there's the additional constraint of having to declare all variables together right after the function declaration, together with any needed constants.

Also important is that in Fortran parentheses are required around if statements and logical operators such as ".and." and ".or.", whereas in Visual Basic they are optional.

The functions declarations and endings are also very different in both programming languages, but they follow a regular pattern in both cases and therefore it's easy to go from one to the other.

And like we did for other languages with built-in min and max functions (thus far, C#, VB.Net, Java, and JavaScript), we are going to make good use of the native min and max functions of Fortran to clean up the code a bit, leaving as usual the code where min and max values are used to evaluate argmin and argmax values.

So here's the Fortran code for the 39 new math functions that were implemented in GeneXproTools 5.0 MR1 "Mini-Release I: What's New?":

    real(8) function gepRamp1(x)
        real(8), intent(in) :: x
        if (x > 0.0) then
            gepRamp1 = x
        else
            gepRamp1 = 0.0
        end if
    end function gepRamp1

    real(8) function gepRamp2(x)
        real(8), intent(in) :: x
        if (x > 0.0) then
            gepRamp2 = 0.0
        else
            gepRamp2 = x
        end if
    end function gepRamp2

    real(8) function gepRamp3(x)
        real(8), intent(in) :: x
        if (x > 0.0) then
            gepRamp3 = 0.0
        else
            gepRamp3 = -x
        end if
    end function gepRamp3

    real(8) function gepRamp4(x)
        real(8), intent(in) :: x
        if (x > 0.0) then
            gepRamp4 = -x
        else
            gepRamp4 = 0.0
        end if
    end function gepRamp4

    real(8) function gepStep1(x)
        real(8), intent(in) :: x
        if (x > 0.0) then
            gepStep1 = 1.0
        else
            gepStep1 = -1.0
        end if
    end function gepStep1

    real(8) function gepStep2(x)
        real(8), intent(in) :: x
        if (x > 0.0) then
            gepStep2 = 1.0
        else
            gepStep2 = 0.0
        end if
    end function gepStep2

    real(8) function gepStep3(x)
        real(8), intent(in) :: x
        if (x >= 1.0) then
            gepStep3 = 1.0
        else if (x <= -1.0) then
            gepStep3 = -1.0
        else
            gepStep3 = x
        end if
    end function gepStep3

    real(8) function gepStep4(x)
        real(8), intent(in) :: x
        if (x >= 1.0) then
            gepStep4 = 1.0
        else if (x <= 0.0) then
            gepStep4 = 0.0
        else
            gepStep4 = x
        end if
    end function gepStep4

    real(8) function gepCL2A(x, y)
        real(8), intent(in) :: x, y
        if ((x > 0.0) .and. (y > 0.0)) then
            gepCL2A = 1.0
        else
            gepCL2A = -1.0
        end if
    end function gepCL2A

    real(8) function gepCL2B(x, y)
        real(8), intent(in) :: x, y
        if ((x >= 0.0) .and. (y < 0.0)) then
            gepCL2B = -1.0
        else
            gepCL2B = 1.0
        end if
    end function gepCL2B

    real(8) function gepCL2C(x, y)
        real(8), intent(in) :: x, y
        if ((x > 1.0) .and. (y < -1.0)) then
            gepCL2C = -1.0
        else
            gepCL2C = 1.0
        end if
    end function gepCL2C

    real(8) function gepCL2D(x, y)
        real(8), intent(in) :: x, y
        if ((x > 0.0) .and. (y > 0.0)) then
            gepCL2D = 1.0
        else
            gepCL2D = 0.0
        end if
    end function gepCL2D

    real(8) function gepCL2E(x, y)
        real(8), intent(in) :: x, y
        if ((x >= 0.0) .and. (y <= 0.0)) then
            gepCL2E = 0.0
        else
            gepCL2E = 1.0
        end if
    end function gepCL2E

    real(8) function gepCL2F(x, y)
        real(8), intent(in) :: x, y
        if ((x > 1.0) .and. (y < -1.0)) then
            gepCL2F = 0.0
        else
            gepCL2F = 1.0
        end if
    end function gepCL2F

    real(8) function gepCL3A(x, y)
        real(8), intent(in) :: x, y
        if ((x > 0.0) .and. (y < 0.0)) then
            gepCL3A = 1.0
        else if ((x < 0.0) .and. (y > 0.0)) then
            gepCL3A = -1.0
        else
            gepCL3A = 0.0
        end if
    end function gepCL3A

    real(8) function gepCL3B(x, y)
        real(8), intent(in) :: x, y
        if ((x >= 1.0) .and. (y >= 1.0)) then
            gepCL3B = 1.0
        else if ((x <= -1.0) .and. (y <= -1.0)) then
            gepCL3B = -1.0
        else
            gepCL3B = 0.0
        end if
    end function gepCL3B

    real(8) function gepCL3C(x, y)
        real(8), intent(in) :: x, y
        if ((x > 0.0) .and. (y > 0.0)) then
            gepCL3C = 1.0
        else if ((x < 0.0) .and. (y < 0.0)) then
            gepCL3C = -1.0
        else
            gepCL3C = 0.0
        end if
    end function gepCL3C

    real(8) function gepMap3A(x, y)
        real(8), intent(in) :: x, y
        real(8), parameter :: SLACK = 10.0
        if (y < (x - SLACK)) then
            gepMap3A = -1.0
        else if (y > (x + SLACK)) then
            gepMap3A = 1.0
        else
            gepMap3A = 0.0
        end if
    end function gepMap3A

    real(8) function gepMap3B(x, y, z)
        real(8), intent(in) :: x, y, z
        real(8) :: minValue
        real(8) :: maxValue
        minValue = min(x,y)
        maxValue = max(x,y)
        
        if (z < minValue) then
            gepMap3B = -1.0
        else if (z > maxValue) then
            gepMap3B = 1.0
        else
            gepMap3B = 0.0
        end if
    end function gepMap3B

    real(8) function gepMap3C(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: minValue
        real(8) :: maxValue
        minValue = min(a,b,c)
        maxValue = max(a,b,c)

        if (d < minValue) then
            gepMap3C = -1.0
        else if (d > maxValue) then
            gepMap3C = 1.0
        else
            gepMap3C = 0.0
        end if
    end function gepMap3C

    real(8) function gepMap4A(x, y)
        real(8), intent(in) :: x, y
        real(8), parameter :: SLACK = 10.0
        if (y < (x - SLACK)) then
            gepMap4A = 0.0
        else if ((y >= (x - SLACK)) .and. (y < x)) then
            gepMap4A = 1.0
        else if ((y >= x) .and. (y < (x + SLACK))) then
            gepMap4A = 2.0
        else if (y >= (x + SLACK)) then
            gepMap4A = 3.0
        end if
    end function gepMap4A

    real(8) function gepMap4B(x, y, z)
        real(8), intent(in) :: x, y, z
        ! evaluate minValue(x,y), maxValue(x,y) and midrange
        real(8) :: minValue
        real(8) :: maxValue
        real(8) :: midrange
        minValue = min(x,y)
        maxValue = max(x,y)
        midrange = (maxValue + minValue)/2.0
        
        if (z < minValue) then
            gepMap4B = 0.0
        else if ((z >= minValue) .and. (z < midrange)) then
            gepMap4B = 1.0
        else if ((z >= midrange) .and. (z < maxValue)) then
            gepMap4B = 2.0
        else if (z >= maxValue) then
            gepMap4B = 3.0
        end if
    end function gepMap4B

    real(8) function gepMap4C(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: minValue, maxValue, midleValue
        integer :: argMin, argMax
        ! evaluate minValue(a,b,c), maxValue(a,b,c) and midleValue(a,b,c)
        !
        ! evaluate minValue(a,b,c)
        minValue = a
        argMin = 0
        if (minValue > b) then
            minValue = b
            argMin = 1
        end if
        if (minValue > c) then
            minValue = c
            argMin = 2
        end if
        ! evaluate maxValue(a,b,c)
        maxValue = a
        argMax = 0
        if (maxValue < b) then
            maxValue = b
            argMax = 1
        end if
        if (maxValue < c) then
            maxValue = c
            argMax = 2
        end if
        ! evaluate midleValue(a,b,c)
        midleValue = c
        if ((0 /= argMin) .and. (0 /= argMax)) then
            midleValue = a
        end if
        if ((1 /= argMin) .and. (1 /= argMax)) then
            midleValue = b
        end if

        if (d < minValue) then
            gepMap4C = 0.0
        else if ((d >= minValue) .and. (d < midleValue)) then
            gepMap4C = 1.0
        else if ((d >= midleValue) .and. (d < maxValue)) then
            gepMap4C = 2.0
        else if (d >= maxValue) then
            gepMap4C = 3.0
        end if
    end function gepMap4C

    real(8) function gepMap5A(x, y)
        real(8), intent(in) :: x, y
        real(8), parameter :: SLACK = 15.0
        if (y < (x - SLACK)) then
            gepMap5A = 0.0
        else if ((y >= (x - SLACK)) .and. (y < (x - SLACK/3.0))) then
            gepMap5A = 1.0
        else if ((y >= (x - SLACK/3.0)) .and. (y < (x + SLACK/3.0))) then
            gepMap5A = 2.0
        else if ((y >= (x + SLACK/3.0)) .and. (y < (x + SLACK))) then
            gepMap5A = 3.0
        else if (y >= (x + SLACK)) then
            gepMap5A = 4.0
        end if
    end function gepMap5A

    real(8) function gepMap5B(x, y, z)
        real(8), intent(in) :: x, y, z
        real(8) :: minValue, maxValue, intervalLength, midpoint1, midpoint2
        ! evaluate minValue(x,y), maxValue(x,y), midpoint1, midpoint2
        minValue = min(x,y)
        maxValue = max(x,y)
        intervalLength = (maxValue - minValue)/3.0
        midpoint1 = minValue + intervalLength
        midpoint2 = minValue + 2.0*intervalLength
        
        if (z < minValue) then
            gepMap5B = 0.0
        else if ((z >= minValue) .and. (z < midpoint1)) then
            gepMap5B = 1.0
        else if ((z >= midpoint1) .and. (z < midpoint2)) then
            gepMap5B = 2.0
        else if ((z >= midpoint2) .and. (z < maxValue)) then
            gepMap5B = 3.0
        else if (z >= maxValue) then
            gepMap5B = 4.0
        end if
    end function gepMap5B

    real(8) function gepMap5C(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: minValue, maxValue, midleValue, midrange1, midrange2
        integer :: argMin, argMax
        ! evaluate minValue(a,b,c), maxValue(a,b,c), midleValue(a,b,c), midrange1, midrange2
        !
        ! evaluate minValue(a,b,c)
        minValue = a
        argMin = 0
        if (minValue > b) then
            minValue = b
            argMin = 1
        end if
        if (minValue > c) then
            minValue = c
            argMin = 2
        end if
        ! evaluate maxValue(a,b,c)
        maxValue = a
        argMax = 0
        if (maxValue < b) then
            maxValue = b
            argMax = 1
        end if
        if (maxValue < c) then
            maxValue = c
            argMax = 2
        end if
        ! evaluate midleValue(a,b,c)
        midleValue = c
        if ((0 /= argMin) .and. (0 /= argMax)) then
            midleValue = a
        end if
        if ((1 /= argMin) .and. (1 /= argMax)) then
            midleValue = b
        end if
        midrange1 = (minValue + midleValue)/2.0
        midrange2 = (midleValue + maxValue)/2.0

        if (d < minValue) then
            gepMap5C = 0.0
        else if ((d >= minValue) .and. (d < midrange1)) then
            gepMap5C = 1.0
        else if ((d >= midrange1) .and. (d < midrange2)) then
            gepMap5C = 2.0
        else if ((d >= midrange2) .and. (d < maxValue)) then
            gepMap5C = 3.0
        else if (d >= maxValue) then
            gepMap5C = 4.0
        end if
    end function gepMap5C

    real(8) function gepMap6A(x, y)
        real(8), intent(in) :: x, y
        real(8), parameter :: SLACK = 10.0
        if (y < (x - SLACK)) then
            gepMap6A = 0.0
        else if ((y >= (x - SLACK)) .and. (y < (x - SLACK/2.0))) then
            gepMap6A = 1.0
        else if ((y >= (x - SLACK/2.0)) .and. (y < x)) then
            gepMap6A = 2.0
        else if ((y >= x) .and. (y < (x + SLACK/2.0))) then
            gepMap6A = 3.0
        else if ((y >= (x + SLACK/2.0)) .and. (y < (x + SLACK))) then
            gepMap6A = 4.0
        else if (y >= (x + SLACK)) then
            gepMap6A = 5.0
        end if
    end function gepMap6A

    real(8) function gepMap6B(x, y, z)
        real(8), intent(in) :: x, y, z
        ! evaluate minValue(x,y), maxValue(x,y), midrange, midpoint1, midpoint2
        real(8) :: minValue, maxValue, midrange, midpoint1, midpoint2
        minValue = min(x,y)
        maxValue = max(x,y)
        midrange = (minValue + maxValue)/2.0
        midpoint1 = (minValue + midrange)/2.0
        midpoint2 = (midrange + maxValue)/2.0
        
        if (z < minValue) then
            gepMap6B = 0.0
        else if ((z >= minValue) .and. (z < midpoint1)) then
            gepMap6B = 1.0
        else if ((z >= midpoint1) .and. (z < midrange)) then
            gepMap6B = 2.0
        else if ((z >= midrange) .and. (z < midpoint2)) then
            gepMap6B = 3.0
        else if ((z >= midpoint2) .and. (z < maxValue)) then
            gepMap6B = 4.0
        else if (z >= maxValue) then
            gepMap6B = 5.0
        end if
    end function gepMap6B

    real(8) function gepMap6C(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: minValue, maxValue, midleValue, midrange1, midrange2
        integer :: argMin, argMax
        ! evaluate minValue(a,b,c), maxValue(a,b,c), midleValue(a,b,c), midrange1, midrange2
        !
        ! evaluate minValue(a,b,c)
        minValue = a
        argMin = 0
        if (minValue > b) then
            minValue = b
            argMin = 1
        end if
        if (minValue > c) then
            minValue = c
            argMin = 2
        end if
        ! evaluate maxValue(a,b,c)
        maxValue = a
        argMax = 0
        if (maxValue < b) then
            maxValue = b
            argMax = 1
        end if
        if (maxValue < c) then
            maxValue = c
            argMax = 2
        end if
        ! evaluate midleValue(a,b,c)
        midleValue = c
        if ((0 /= argMin) .and. (0 /= argMax)) then
            midleValue = a
        end if
        if ((1 /= argMin) .and. (1 /= argMax)) then
            midleValue = b
        end if
        ! evaluate midrange1 and midrange2
        midrange1 = (minValue + midleValue)/2.0
        midrange2 = (midleValue + maxValue)/2.0

        if (d < minValue) then
            gepMap6C = 0.0
        else if ((d >= minValue) .and. (d < midrange1)) then
            gepMap6C = 1.0
        else if ((d >= midrange1) .and. (d < midleValue)) then
            gepMap6C = 2.0
        else if ((d >= midleValue) .and. (d < midrange2)) then
            gepMap6C = 3.0
        else if ((d >= midrange2) .and. (d < maxValue)) then
            gepMap6C = 4.0
        else if (d >= maxValue) then
            gepMap6C = 5.0
        end if
    end function gepMap6C

    real(8) function gepECL3A(x, y, z)
        real(8), intent(in) :: x, y, z
        if ((y > x) .and. (z < x)) then
            gepECL3A = 1.0
        else if ((y < x) .and. (z > x)) then
            gepECL3A = -1.0
        else
            gepECL3A = 0.0
        end if
    end function gepECL3A

    real(8) function gepECL3B(x, y, z)
        real(8), intent(in) :: x, y, z
        if ((y > x) .and. (z > x)) then
            gepECL3B = 1.0
        else if ((y < x) .and. (z < x)) then
            gepECL3B = -1.0
        else
            gepECL3B = 0.0
        end if
    end function gepECL3B

    real(8) function gepECL3C(x, y, z)
        real(8), intent(in) :: x, y, z
        if ((y >= x) .and. (z >= x)) then
            gepECL3C = 1.0
        else if ((y <= -x) .and. (z <= -x)) then
            gepECL3C = -1.0
        else
            gepECL3C = 0.0
        end if
    end function gepECL3C

    real(8) function gepECL3D(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: minValue, maxValue
        minValue = min(a,b)
        maxValue = max(a,b)

        if ((c >= maxValue) .and. (d >= maxValue)) then
            gepECL3D = 1.0
        else if ((c <= minValue) .and. (d <= minValue)) then
            gepECL3D = -1.0
        else
            gepECL3D = 0.0
        end if
    end function gepECL3D

    real(8) function gepAMin2(x, y)
        real(8), intent(in) :: x, y
        if (x < y) then
            gepAMin2 = 0.0
        else
            gepAMin2 = 1.0
        end if
    end function gepAMin2

    real(8) function gepAMin3(x, y, z)
        real(8), intent(in) :: x, y, z
        real(8) :: temp
        real(8) :: argMin
        temp = x
        argMin = 0.0    
        if (temp >= y) then
            temp = y
            argMin = 1.0
        end if
        if (temp >= z) then
            argMin = 2.0
        end if
        gepAMin3 = argMin
    end function gepAMin3

    real(8) function gepAMin4(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: temp
        real(8) :: argMin
        temp = a
        argMin = 0.0    
        if (temp >= b) then
            temp = b
            argMin = 1.0
        end if
        if (temp >= c) then
            temp = c
            argMin = 2.0
        end if
        if (temp >= d) then
            argMin = 3.0
        end if
        gepAMin4 = argMin
    end function gepAMin4

    real(8) function gepAMax2(x, y)
        real(8), intent(in) :: x, y
        if (x >= y) then
            gepAMax2 = 0.0
        else
            gepAMax2 = 1.0
        end if
    end function gepAMax2

    real(8) function gepAMax3(x, y, z)
        real(8), intent(in) :: x, y, z
        real(8) :: temp
        real(8) :: argMax
        temp = x
        argMax = 0.0    
        if (temp < y) then
            temp = y
            argMax = 1.0
        end if
        if (temp < z) then
            argMax = 2.0
        end if
        gepAMax3 = argMax
    end function gepAMax3

    real(8) function gepAMax4(a, b, c, d)
        real(8), intent(in) :: a, b, c, d
        real(8) :: temp
        real(8) :: argMax
        temp = a
        argMax = 0.0    
        if (temp < b) then
            temp = b
            argMax = 1.0
        end if
        if (temp < c) then
            temp = c
            argMax = 2.0
        end if
        if (temp < d) then
            argMax = 3.0
        end if
        gepAMax4 = argMax
    end function gepAMax4

Author

Comments

There are currently no comments on this article.

Comment

your_ip_is_blacklisted_by sbl.spamhaus.org

← Older Newer →