The SIM2D Procedure |
The subregion to be considered is the southeast corner of the field, which is a square region with a length of 40 distance units (in thousands of feet). First, you input the thickness data as the following DATA step shows:
data thick; input East North Thick @@; label Thick='Coal Seam Thickness'; datalines; 0.7 59.6 34.1 2.1 82.7 42.2 4.7 75.1 39.5 4.8 52.8 34.3 5.9 67.1 37.0 6.0 35.7 35.9 6.4 33.7 36.4 7.0 46.7 34.6 8.2 40.1 35.4 13.3 0.6 44.7 13.3 68.2 37.8 13.4 31.3 37.8 17.8 6.9 43.9 20.1 66.3 37.7 22.7 87.6 42.8 23.0 93.9 43.6 24.3 73.0 39.3 24.8 15.1 42.3 24.8 26.3 39.7 26.4 58.0 36.9 26.9 65.0 37.8 27.7 83.3 41.8 27.9 90.8 43.3 29.1 47.9 36.7 29.5 89.4 43.0 30.1 6.1 43.6 30.8 12.1 42.8 32.7 40.2 37.5 34.8 8.1 43.3 35.3 32.0 38.8 37.0 70.3 39.2 38.2 77.9 40.7 38.9 23.3 40.5 39.4 82.5 41.4 43.0 4.7 43.3 43.7 7.6 43.1 46.4 84.1 41.5 46.7 10.6 42.6 49.9 22.1 40.7 51.0 88.8 42.0 52.8 68.9 39.3 52.9 32.7 39.2 55.5 92.9 42.2 56.0 1.6 42.7 60.6 75.2 40.1 62.1 26.6 40.1 63.0 12.7 41.8 69.0 75.6 40.1 70.5 83.7 40.9 70.9 11.0 41.7 71.5 29.5 39.8 78.1 45.5 38.7 78.2 9.1 41.7 78.4 20.0 40.8 80.5 55.9 38.7 81.1 51.0 38.6 83.8 7.9 41.6 84.5 11.0 41.5 85.2 67.3 39.4 85.5 73.0 39.8 86.7 70.4 39.6 87.2 55.7 38.8 88.1 0.0 41.6 88.4 12.1 41.3 88.4 99.6 41.2 88.8 82.9 40.5 88.9 6.2 41.5 90.6 7.0 41.5 90.7 49.6 38.9 91.5 55.4 39.0 92.9 46.8 39.1 93.4 70.9 39.7 55.8 50.5 38.1 96.2 84.3 40.3 98.2 58.2 39.5 ;
PROC SIM2D is run on the entire data set for conditioning, while the simulation grid covers only this subregion. It is convenient to be able to vary the seed, the grid increment, and the number of simulations performed. The following macro implements the computation of the percent area exceeding the cutoff value by using the seed, the grid increment, and the number of simulated realizations as macro arguments.
Within the macro, the data set produced by PROC SIM2D is transposed with PROC TRANSPOSE so that each grid location is a separate variable. The MEANS procedure is then used to average the simulated value at each grid point over all realizations. It is this average that is compared to the cutoff value. The last DATA step does the comparison, uses an appropriate loop to determine the percent of the grid locations that exceed this cutoff value, and writes the results to the listing file in the form of a report. This sequence is implemented with the following statements:
/* Construct macro for conditional simulation ------------------*/ %let cc0=7.2881; %let aa0=30.6239; %let ngt=1e-8; %let form=gauss; %let cut=39.7; %macro area_sim(seed=,nr=,ginc=); %let ngrid=%eval(40/&ginc+1); %let tgrid=%eval(&ngrid*&ngrid); proc sim2d data=thick outsim=sim1; coordinates xc=east yc=north; simulate var=thick numreal=&nr seed=&seed scale=&cc0 range=&aa0 nugget=&ngt form=&form; mean 40.1173; grid x=60 to 100 by &ginc y= 0 to 40 by &ginc; run; proc transpose data=sim1 out=sim2 prefix=sims; by _iter_; var svalue; run; proc means data=sim2 noprint n mean; var sims1-sims&tgrid; output out=msim n=numsim mean=ms1-ms&tgrid; run; data _null_; file print; array simss ms1-ms&tgrid; set msim; cflag=0; do ss=1 to &tgrid; tempv=simss[ss]; if simss[ss] > &cut then do; cflag + 1; end; end; area_per=100*(cflag/&tgrid); put // +5 'Conditional Simulation of Coal Seam' ' Thickness for Subregion'; put / +5 'Subregion is South-East Corner 40 by 40 distance units'; put / +5 "Seed:&seed" +2 "Grid Increment:&ginc"; put / +5 "Total Number of Grid Points:&tgrid" +2 "Number of Simulations:&nr"; put / +5 "Percent of Subregion Exceeding Cutoff of %left(&cut) ft.:" +2 area_per 5.2; run; %mend area_sim;
In the following statement, you invoke the macro three times. Each time, the macro is invoked with a different seed and combination of the grid increment and number of simulations. The macro is first invoked with a relatively coarse grid (grid increment of 10 distance units) and a small number of realizations (5). The output of this conditional simulation is shown in Output 79.1.1.
/* Execute macro for coarse grid -------------------------------*/ %area_sim(seed=12345,nr=5,ginc=10);
The next invocation, in the following statement, uses a finer grid and 50 realizations. The output of the second conditional simulation is shown in Output 79.1.2.
/* Execute macro for fine grid and fewer simulations -----------*/ %area_sim(seed=54321,nr=50,ginc=1);
The final invocation, in the following statement, uses the same grid increment and 500 realizations. The output of this conditional simulation is shown in Output 79.1.3.
/* Execute macro for fine grid and more simulations ------------*/ %area_sim(seed=655311,nr=500,ginc=1);
The results from the preceding simulations indicate that 76% to 77% of the subregion exceeds the cutoff value.
Note also that the number of grid points in the simulation increases with the square of the decrease in the grid increment, leading to long CPU processing times. Increasing the number of realizations results in a linear increase in processing times. Hence, using as coarse a grid as possible allows for more realizations and experimentation with different seeds.
Copyright © 2009 by SAS Institute Inc., Cary, NC, USA. All rights reserved.