subroutine ph_eqfac (npt,method,k,hp,alk,hco3,hs,boh4,fco2,fco3,
     &             fhco3,fh2s,fhs,fboh3,fboh4,fahco3,faco3,faboh4,fahs)
C-------------------------------------------------------------------
C PH_EQFAC $Revision: 1.0 $ $Date: 1998/12/07 $
C          Copyright (C) GEOMAR, Roger Luff 1998.
C
C USAGE: 
C     INTEGER npt
C     parameter (npt=10)
C     double precision fco2(npt),fco3(npt),fhco3(npt)
C     double precision fh2s(npt),fhs(npt),fboh3(npt),fboh4(npt)
C     double precision fahco3(npt),faco3(npt),faboh4(npt),fahs(npt)
C     call ph_eqfac (npt,1,k,hp,alk,hco3,hs,boh4,fco2,fco3,fhco3,
C    &       fh2s,fhs,fboh3,fboh4,fahco3,faco3,faboh4,fahs)
C   SEE ph_test.f FOR DETAILS OF USING THIS PROGRAM!!!!!!!!!!!!!!!
C
C DESCRIPTION:
C   This subroutine calculates the mol fractions of the individual
C   species with respect to the sun concentrations alkalinity, TCO2
C   TH2S, TBOH for the description of the transport processes in a
C   numerical model (see reference for details)
C
C INPUT:
C   npt        = number of layers in the model
C   method     = type of equilibrium calculation: 
C                1 alkalinity conservation approach 
C                2 advancement approach
C                3 charge conservation approach
C   k          = equilibrium constants from ph_thermo
C   hp         = H+ concentration [mmol/l]
C   hco3       = Equilibrium concentration of HCO3 [mmol/l]
C   hs         = Equilibrium concentration of HS [mmol/l]
C
C OUTPUT:
C   boh4       = Equilibrium concentration of BOH4 [mmol/l]
C   fco2       = mole fraction of co2 the respect to TCO2
C   fco3       = mole fraction of co3 the respect to TCO2
C   fhco3      = mole fraction of hco3 the respect to TCO2
C   fh2s       = mole fraction of h2s the respect to TH2S
C   fhs        = mole fraction of hs the respect to TH2S
C   fboh3      = mole fraction of boh3 the respect to TBOH
C   fboh4      = mole fraction of boh4 the respect to TBOH
C   fahco3     = mole fraction of hco3 the respect to alkalinity
C   faco3      = mole fraction of co3 the respect to alkalinity
C   faboh4     = mole fraction of boh4 the respect to alkalinity
C   fahs        = mole fraction of hs the respect to alkalinity
C
C AUTHOR:
C   Roger Luff and Matthias Haeckel 98-12-07  (rluff@gmx.de)
C   GEOMAR Research Centre for Marine Geosciences
C   Wischhofstr. 1-3, D-24148 Kiel, Germany
C   visumod@technologist.com
C
C DISCLAIMER:
C   This software is provided "as it is" without warranty of any kind.
C
C REFERENCES:
C   Luff, Haeckel and Wallmann (submitted) "Robust and fast FORTRAN
C   and MATLAB libraries to calculate pH distributions in
C   marine systems", Computers & Geosciences.
C-------------------------------------------------------------------

      implicit none
      integer npt,method
      double precision k(12)
      double precision alk(npt),hco3(npt),hs(npt),boh4(npt)
      double precision hp(npt),fco2(npt),fco3(npt),fhco3(npt)
      double precision fh2s(npt),fhs(npt),fboh3(npt),fboh4(npt)
      double precision fahco3(npt),faco3(npt),faboh4(npt),fahs(npt)
      integer i

      do i=1,npt      
        if (method .eq. 1) then
C-------------------------------------------------------------------
c         Factors for alkalinity conservation approach
C-------------------------------------------------------------------
          fhco3(i) = 1.0/(1.0+hp(i)/k(2)+k(3)/hp(i))
          fco2(i) = 1.0/(1.0+k(2)/hp(i)+k(2)*k(3)/(hp(i)**2))
          fco3(i) = 1.0-(fhco3(i)+fco2(i))
          fahco3(i) = hco3(i)/alk(i)
          fahs(i) = hs(i)/alk(i)
          faboh4(i) = boh4(i)/alk(i)
          faco3(i) = 1.0-(fahco3(i)+fahs(i)+faboh4(i))
        endif
        if (method.ne.2) then
C-------------------------------------------------------------------
c         No Factors are necessary for advancement approach
C-------------------------------------------------------------------
          fboh3(i) = (1.0/(1.0+k(4)/hp(i)))
          fboh4(i) = 1.0-fboh3(i)
          fhs(i) = (1.0/(1.0+hp(i)/k(10)))
          fh2s(i) = 1.0-fhs(i)
        endif
        if (method .eq. 3) then
C-------------------------------------------------------------------
c         Factors for charge conservation approach
C-------------------------------------------------------------------
          fhco3(i) = (1.0/(1.0+hp(i)/k(2)+k(3)/hp(i)))
          fco2(i) = (1.0/(1.0+k(2)/hp(i)+k(2)*k(3)/(hp(i)**2)))
          fco3(i) = 1.0-(fhco3(i)+fco2(i))
          fhs(i) = (1.0/(1.0+hp(i)/k(10)))
        endif
      enddo
      return
      end