private double surfacePotentialBisection(double pLow, double pHigh, double sigma, int n){
        double surPot = 0.0D;   // temporary values of the surface potential

        // Bisection method converging on sigma calculated = true sigma
        boolean testC = true;
        int testCiter = 0;
        int testCmax = 10;
        double pDiff = pHigh-pLow;
        double cLow = cFunction(pLow, sigma);
        double cHigh = cFunction(pHigh, sigma);
        while(testC){
            if(pHigh*pLow>0.0D){
                testCiter++;
                if(testCiter>testCmax)throw new IllegalArgumentException("root not bounded after " + testCiter + "expansions");
                pLow  -= pDiff;
                pHigh += pDiff;
                cLow = cFunction(pLow, sigma);
                cHigh = cFunction(pHigh, sigma);
            }
            else{
                testC=false;
            }
        }
        double check = Math.abs(sigma)*1e-6;
        boolean test = true;
        double cMid = 0.0D;
        int nIter = 0;
        while(test){
            surPot = (pLow + pHigh)/2.0D;
            cMid = this.cFunction(surPot, sigma);
            if(Math.abs(cMid)<=check){
                test = false;
            }
            else{
                nIter++;
                if(nIter>10000){
                    System.out.println("Class: GouyChapmanStern\nMethod: surfacePotential"+n +"\nnumber of iterations exceeded in bisection\ncurrent value of sigma returned");
                    test = false;
                }
                if(cMid*cHigh>0){
                    pHigh = surPot;
                    cHigh = cMid;
                }
                else{
                    pLow = surPot;
                    cLow = cMid;
                }
            }
        }

        return surPot;
    }