CUDA-Based Implementation of GSLIB: The Geostatistical Software - PowerPoint PPT Presentation

CUDA-Based Implementation of GSLIB: The Geostatistical Software Library Daniel Baeza Oscar Peredo dabaeza@alges.cl operedo@alges.cl

The Mining Process Exploration Evaluation Planning Operation

GSLIB: The Geostatistical Software Library �

GSLIB: The Geostatistical Software Library • GSLIB is a software package composed by a set of utilities and applications related with geostatistics • Full implemented in Fortran 77/90 • Run in OSX, Linux and Windows • Widely used for academics, researchers, engineers GSLIB: Geostatistical Software Library and user's guide (1998) Deutsch, Clayton V, Journel, André G

Variogram calculation with GSLIB • gamv is the GSLIB variogram calculation method • It’s a fundamental tool in geostatistics • Allow to quantify the spatial variability of a variable • Used in geostatistical estimation and simulation Normal Scores Semivariogram Low Solubil Normal Scores Semivariogram High Solubi 1.20 1.20 • High computational cost 1.00 1.00 .80 .80 γ γ .60 .60 .40 .40 .20 .20 .00 .00 0. 100. 200. 300. 400. 500. 0. 100. 200. 300. 400. 500. Distance Distance

Variogram calculation Z(u) 2 ! (h) = 1 ∑ [z(u) - z(u + h)] 2 || || = h a - b b N(h) z Z(u + h) a y x ! (h) more variability less variability h

Variogram computation � Sequential & Parallel Implementations

Sequential implementation Input : • ( VR , Ω ): sample data values VR ( m columns) defined in a 3D domain of coordinates Ω Setup parameters • nvar : number of variables ( nvar ≤ m ) & • nlag : number of lags • h : lag separation distance Load data • ndir : number of directions • h 1 , . . . , h ndir : directions • τ 1 , . . . , τ ndir : geometrical tolerance parameters • nvarg : number of variograms • ( ivtype 1 , ivtail 1 , ivhead 1 ) , . . . , ( ivtype nvarg , ivtail nvarg , ivhead nvarg ): variogram types 1 Read input parameter file; 2 Read sample data values file; 3 β ← zeros ( nvar × nlag × ndir × nvarg ); 4 for i ∈ { 1 , . . . , | Ω |} do Main computation: Loop over pairs of points for j ∈ { i , . . . , | Ω |} do 5 for id ∈ { 1 , . . . , ndir } do 6 for iv ∈ { 1 , . . . , nvarg } do 7 for il ∈ { 1 , . . . , nlag } do 8 p i = ( x i , y i , z i ) ∈ Ω ; 9 p j = ( x j , y j , z j ) ∈ Ω ; 10 if ( p i , p j ) satisfy tolerances τ id and || p i − p j || ≈ h id × il × h then 11 Save ( VR i , ivhead iv or VR i , ivtail iv ) and ( VR j , ivhead iv or VR j , ivtail iv ) according to variogram ivtype iv into β ; 12 Read computation results 13 γ ← build variogram using statistics β Write γ in the output file Output : Output file with γ values

Sequential implementation Input : • ( VR , Ω ): sample data values VR ( m columns) defined in a 3D domain of coordinates Ω • nvar : number of variables ( nvar ≤ m ) • nlag : number of lags • h : lag separation distance • ndir : number of directions • h 1 , . . . , h ndir : directions • τ 1 , . . . , τ ndir : geometrical tolerance parameters • nvarg : number of variograms • ( ivtype 1 , ivtail 1 , ivhead 1 ) , . . . , ( ivtype nvarg , ivtail nvarg , ivhead nvarg ): variogram types 1 Read input parameter file; 2 Read sample data values file; 3 β ← zeros ( nvar × nlag × ndir × nvarg ); 4 for i ∈ { 1 , . . . , | Ω |} do for j ∈ { i , . . . , | Ω |} do 5 for id ∈ { 1 , . . . , ndir } do 6 for iv ∈ { 1 , . . . , nvarg } do 7 for il ∈ { 1 , . . . , nlag } do 8 p i = ( x i , y i , z i ) ∈ Ω ; 9 p j = ( x j , y j , z j ) ∈ Ω ; 10 if ( p i , p j ) satisfy tolerances τ id and || p i − p j || ≈ h id × il × h then 11 Save ( VR i , ivhead iv or VR i , ivtail iv ) and ( VR j , ivhead iv or VR j , ivtail iv ) according to variogram ivtype iv into β ; 12 13 γ ← build variogram using statistics β Write γ in the output file Output : Output file with γ values

Sequential implementation Input : • ( VR , Ω ): sample data values VR ( m columns) defined in a 3D domain of coordinates Ω • nvar : number of variables ( nvar ≤ m ) • nlag : number of lags • h : lag separation distance • ndir : number of directions • h 1 , . . . , h ndir : directions • τ 1 , . . . , τ ndir : geometrical tolerance parameters • nvarg : number of variograms • ( ivtype 1 , ivtail 1 , ivhead 1 ) , . . . , ( ivtype nvarg , ivtail nvarg , ivhead nvarg ): variogram types 1 Read input parameter file; 2 Read sample data values file; 3 β ← zeros ( nvar × nlag × ndir × nvarg ); 4 for i ∈ { 1 , . . . , | Ω |} do for j ∈ { i , . . . , | Ω |} do 5 for id ∈ { 1 , . . . , ndir } do 6 for iv ∈ { 1 , . . . , nvarg } do 7 for il ∈ { 1 , . . . , nlag } do 8 p i = ( x i , y i , z i ) ∈ Ω ; 9 p j = ( x j , y j , z j ) ∈ Ω ; 10 if ( p i , p j ) satisfy tolerances τ id and || p i − p j || ≈ h id × il × h then 11 Save ( VR i , ivhead iv or VR i , ivtail iv ) and ( VR j , ivhead iv or VR j , ivtail iv ) according to variogram ivtype iv into β ; 12 13 γ ← build variogram using statistics β Write γ in the output file Output : Output file with γ values

• • ( ivtype 1 , ivtail 1 , ivhead 1 ) , . . . , ( ivtype nvarg , ivtail nvarg , ivhead nvarg ): variogram types Sequential implementation 1 Read input parameter file; 2 Read sample data values file; 3 β ← zeros ( nvar × nlag × ndir × nvarg ); 4 for i ∈ { 1 , . . . , | Ω |} do for j ∈ { i , . . . , | Ω |} do 5 for id ∈ { 1 , . . . , ndir } do 6 compute statistics ( p i , p j ) for iv ∈ { 1 , . . . , nvarg } do 7 store result in shared memory for il ∈ { 1 , . . . , nlag } do 8 p i = ( x i , y i , z i ) ∈ Ω ; 9 p j = ( x j , y j , z j ) ∈ Ω ; 10 if ( p i , p j ) satisfy tolerances τ id and || p i − p j || ≈ h id × il × h then 11 Save ( VR i , ivhead iv or VR i , ivtail iv ) and ( VR j , ivhead iv or VR j , ivtail iv ) according to variogram ivtype iv into β ; 12 13 γ ← build variogram using statistics β Write γ in the output file Output : Output file with γ values

• • ( ivtype 1 , ivtail 1 , ivhead 1 ) , . . . , ( ivtype nvarg , ivtail nvarg , ivhead nvarg ): variogram types Sequential implementation 1 Read input parameter file; 2 Read sample data values file; 3 β ← zeros ( nvar × nlag × ndir × nvarg ); 4 for i ∈ { 1 , . . . , | Ω |} do for j ∈ { i , . . . , | Ω |} do 5 for id ∈ { 1 , . . . , ndir } do 6 compute statistics ( p i , p j ) for iv ∈ { 1 , . . . , nvarg } do 7 store result in shared memory for il ∈ { 1 , . . . , nlag } do 8 p i = ( x i , y i , z i ) ∈ Ω ; 9 p j = ( x j , y j , z j ) ∈ Ω ; 10 if ( p i , p j ) satisfy tolerances τ id and || p i − p j || ≈ h id × il × h then 11 Save ( VR i , ivhead iv or VR i , ivtail iv ) and ( VR j , ivhead iv or VR j , ivtail iv ) according to variogram ivtype iv into β ; 12 13 γ ← build variogram using statistics β Write γ in the output file Output : Output file with γ values STEP N STEP 1 STEP 2 STEP 3 ...

Parallel implementation • 1 id x = blockId.x*blockDim.x + threadId.x /* x threads coord in the GPU grid */ GPU Kernel 2 id y = blockId.y*blockDim.y + threadId.y /* y threads coord in the GPU grid */ 3 Set and Initialize shared memory in the block syncthreads() 4 5 iter x = id x 6 iter y = id y 7 while ( iter x & iter y ∈ ChunkPoints ) /* Chunk points belongs thread (x,y) */ 8 do j = iter x + | Ω | 9 2 i = iter y 10 compute statistics ( p i , p j ) 11 store result in shared memory via atomic functions 12 if ( iter x > iter y ) then 13 i = iter y 14 j = iter x 15 Each thread compute compute statistics ( p i , p j ) 16 store result in shared memory via atomic functions 17 else only a couple of correlation values 18 if ( iter x == iter y ) then 19 i = iter y 20 j = iter y 21 compute statistics ( p i , p j ) 22 store result in shared memory via atomic functions 23 i = iter x + | Ω | 24 2 j = iter y + | Ω | 25 2 compute statistics ( p i , p j ) 26 store result in shared memory via atomic functions 27 up date ( iter x , iter y ) 28 syncthreads() 29 30 save Statistics values that are in shared memory into globa memory via atomic functions Output : Array β

Parallel implementation • 1 id x = blockId.x*blockDim.x + threadId.x /* x threads coord in the GPU grid */ 2 id y = blockId.y*blockDim.y + threadId.y /* y threads coord in the GPU grid */ 3 Set and Initialize shared memory in the block syncthreads() 4 5 iter x = id x 6 iter y = id y 7 while ( iter x & iter y ∈ ChunkPoints ) /* Chunk points belongs thread (x,y) */ 8 do j = iter x + | Ω | 9 2 i = iter y 10 compute statistics ( p i , p j ) 11 store result in shared memory via atomic functions 12 if ( iter x > iter y ) then 13 i = iter y 14 j = iter x 15 compute statistics ( p i , p j ) 16 store result in shared memory via atomic functions 17 else 18 if ( iter x == iter y ) then 19 i = iter y 20 j = iter y 21 compute statistics ( p i , p j ) 22 store result in shared memory via atomic functions 23 i = iter x + | Ω | 24 2 j = iter y + | Ω | 25 2 compute statistics ( p i , p j ) 26 store result in shared memory via atomic functions 27 up date ( iter x , iter y ) 28 syncthreads() 29 30 save Statistics values that are in shared memory into globa memory via atomic functions Output : Array β