Figure 2. Semilog plot of the function x = erfc (y) showing the transformed variable (The discontinuity which occurs at z

z.

=

2, and the fact that y= 3.76656 for 0 < z

< 2 reflect the value of y calculated by the INERF subroutine rather than the true value of y as explained in sec. 3.4.)

Figure 3.

Calculation of y

=

erfc-1 (x) by means of a parabolic fit. (The fit is made in terms of the variable z which is calculated from the given x by eq (17).)

Thus knowing a, b, and c, y can be calculated from z, the transformed value of x, by

(24)

(25)

Two forms of data input are employed. The C-V data are read from a file which is titled CVIN. The file CVIN is entered into the computer before CV1 is executed. Each line of CVIN consists of a line number, the capacitance in picofarads, and the applied voltage in volts, where reverse bias is considered positive. Data are listed in order of decreasing capacitance. A sample CVIN file, based on data reported in [10], is given in Appendix C. The second form of data input involves entry of scaling and diffused layer parameters as called for by INPUT statements. The program statements which request these data are lines 2514 to 2585 (see listing in Appendix A). The relevant parameters are (1) G1 or 1ST LOG, the minimum value of dopant density which is within the range of the plot (that is, the dopant density at the origin of the plot), (2) G2 or CHAR/DEC, the number of teletypewriter carriage spaces between decades of dopant density (the plot is semi-log), (3) G3 or LINES/UM, the number of teletypewriter line feeds to represent 1 um on the plot (all positive values are permitted), (4) D or DIAMETER, the diameter of the diode in mils, (5) X5 or JCT DEPTH, the junction depth of the diffused layer in micrometres, (6) NO or SURFACE DENS, the density of diffusant atoms at the surface of the diffused layer (as measured for instance by junction depth and sheet resistivity on a pilot slice [15]) in (centimetres)-3, and (7) N3 or ESTIMATED NB, the background dopant density within the diffused layer in (centimetres)-3. An example of such data input is given in the sample calculation of Appendix C.

3.2. Data Output

Data output is in the form of a plot of dopant density versus distance from the junction, N(B) vs. B. Following this plot, final values of A and N(A) representing the maximum penetration of the depletion region into the diffused layer and the net dopant density at that point, are printed out (lines 3990 to 3992) along with the percent peripheral capacitance correction at the last data point (i.e., the maximum peripheral capacitance correction), and an indication of whether the range of the INERF subroutine has been exceeded (see sec. 3.4.). Finally a true Gaussian profile is calculated and plotted for comparison with the experimental plot. The experimental N(B) vs. B values are written into a file named NBVSB, and the numerical values of the true Gaussian plot are written into a file named IDEAL. These can be listed if desired after the program has run as is done in Appendix C.

3.3. The ERF Subroutine

The ERF subroutine consists of lines 1400 to 1565 in the listing of Appendix A. When this subroutine is called it operates on the variable E, calculating El = erf (E) and J1 = erfc (E) where E corresponds to x in eqs (13) through (15). If E < 10−32, El and J1 are set equal to 0 and 1 respectively (lines 1400 to 1406). If E > 9, El and J1 are set to equal 1 and 0 respectively (lines 1410 to 1416). If 10-3 -32 < E < 1, El and J1 are calculated from eqs (13) and (14). In lines 1430 to 1434 the first 10 terms (to N = 9) of the power series are summed. Additional terms would contribute to the seventh and higher significant figures which is beyond the range of BASIC. In line 1440 the series is multiplied by the pre-factor indicated in eq (13) to get El = erf (E).

If 9 > E > 1, J1 = erfc (E) is calculated from eq (15). The first approximation to the continued fraction is of the form 2x/(2x2 + 1), as calculated in lines 1500 to 1508 and represented by F5. A second approximation to the continued fraction includes the term 2x2 + 5 in eq (15) and is calculated in lines 1510, 1512, and 1520 to 1526 and represented by F9. The first and second approximations are compared in lines 1528 and 1530. If the first and second approximations differ by 1 part in 106 or more, a third approximation is calculated in lines 1538 to 1546 and 1520 to 1526 to include the term 2x2 + 9 in eq (15). The second and third approximations are then compared to determine whether the continued fraction has converged.

Continued iterations are made, each time comparing the present with the previous approximation, until a relative difference of less than 1 part in 106 is obtained. In line 1550 the continued fraction is multiplied by the pre-factor indicated in eq (15) to get J1 = erfc (E). Each approximation of the continued fraction is expressed as a quotient of two expressions (line 1526). As the order of the approximation increases, the number of terms which are multiplied together to construct the two expressions increases. The two expressions could become excessively large and cause overflow. This is prevented by dividing several terms by E3 = 1025 whenever the denominator of this quotient exceeds 1025. This is done in lines 1532 to 1537. The division by this scale factor does not affect the result since it is applied to both numerator and denominator of the quotient.

3.4. The INEFF Subroutine

The INERF subroutine consists of lines 1000 to 1080 in the listing of Appendix A. When this subroutine is called, it operates on the variable H9 and calculates E = erfc1 (H9) by a table look-up approach. Note that H9 in the program corresponds to x in eq (16) and E corresponds to y.

If H9 > 1 (x > 1 in figure 2) the inverse complementary error function is not defined; E is set equal to zero and the flag Z5 (which is printed when the main program is run) is set equal to 2 (lines 1000 to 1003). If H9 = 1, E is set equal to zero (lines 1004 to 1007). The transformed variable HO which corresponds to z is calculated in line 1010 by eq (17). If HO is less than 2 (corresponding to H9 less than 1.39 × 10-7) E is set equal to 3.76656 and flag 25 is set equal to 1 in lines 1012 to 1018. This lower limit arises because (1) the transformation of eq (17) was defined to be applicable for x in the range 10-7 to 1, (2) the value of x corresponding to z = O does not appear in the data table which runs form z = 1 to 100, and (3) the parabolic fit requires that z not be less than z2 (eq (18)) so that interpolation cannot be performed for values of x (H9) for which z (HO) is less than 2. Factors (2) and (3) therefore serve to make the lower limit of H9 be 1.39 × 10-7 rather than 10-7.

In lines 1020 to 1045, Z1, Z2, Z3, Y1, y2, and y3 are defined as H1, H2, H3, H5, H6, and H7, respectively (see figure 3). The coefficients of the parabolic fit a, b, and c or Il, 13, and 14 are calculated by means of eqs (23), (24), and (22), respectively in lines 1050 to 1065. The unknown E (y) is then calculated in line 1070 as in eq (25).

3.5. The TABLE Data File

The TABLE data file consists of lines 150 to 188 in the listing of Appendix A. As discussed in sections 2.4. and 3.4. and referring to eqs (16) and (17), TABLE is a list of y (E) values which correspond to integer values of z (HO) from 1 to 100 inclusive and which facilitate the computation of the function E = erfc-1 (H9) for H9 values from 1.39 × 10 7 to 1.

The pur

The PLOT subroutine consists of lines 2150 to 2460 in the listing of Appendix A. pose of the PLOT subroutine is to produce a semilog plot of dopant density N(B) across the paper of the teletypewriter printout versus distance from the junction B down the paper. The subroutine is also used to plot the true Gaussian profile for comparison with the experimental profile.

The PLOT subroutine can be understood with the help of the flow chart of figure 4 and the partial table of symbols in table 1. Some of the line numbers in the main program are also found in table 1. The flags Z3, Z6, and Z7 are all initially zero. The FLOT subroutine operates on the quantities X1 and Y1 which are calculated from B and N2 (where N2 corresponds to N(B)) in lines 3240 and 3270. These in turn are derived from the C-V data (lines 3000 to 3230) in the first plot and calculated values (lines 4010 to 4050) in the second plot. The first X1, Yl pair (Z6 = 0) is stored as X2, Y2. The Y1 values of succeeding X1, Y1 pairs are accumulated in Y2 until an X1, Y1 pair is input for which X2 > X1, that is, which should be plotted at a farther position along the depth axis. When this occurs Y2 is averaged and integerized. If no data points have yet been plotted (X0 = 0), + signs are printed at decade intervals along the dopant density axis and XO is set equal to 1 before proceeding further. The depth X2 is compared with the line space on the plot XO. If X2 > XO, a series of tick marks in the form of + signs are printed along the depth axis until XO is advanced to X2 (where the data point is to be plotted). The origin of the depth axis is the first of these