GEOL. CROAT. 49/2 145 - 1531 12 Figs. 1 
1 

ZAGREB 1996 

System Architecture for 3D Gravity Modelling 

Franjo SUMANOVAC ', [po L. RITSEMA' and Bert H. BRIL' 

Key word s: Gravity modelling, Com plex models, Soft­
ware archi tecture. 

Abstract 
A fl exi ble so ft ware architecture fo r grav ity modell ing is cS I:lb­

lished and the ad vantage is discussed of hav ing seve ral alternative 
programs 10 handle complex 3D models. The flexible architecture 
consists of four pa ris. implcmclUcd in a d istributed computer environ­
ment : (ile three-dimensional model bui lders and v isual izers (GOCAO 
software. version 7.0), the model representation translators (COCA O 
software or GEOMOD sonwarc), the forward simulation algorit hms 
of gravimetric data ('lpp lyi ng Talwan i-Ewing and GtHze- Luhmcye r 
me thods in the finite-c lement representation class), and thc invcrsion 
(mode l upd;lt ing) sc heme manager based on the Co rde ll - li enderso n 
in ve rsio n proccdure. 

A good so ftware architec ture shou ld at least keep the modcl 
build ing and updating so rt ware separate from the forward simu lation 
software. Inversion schemes can then bc real ized by commun icat ion 
be tween thc two part s of software. 

Several synthetic cascs arc shown 10 demonstrate the lise and the 
capabi lity of the architecture and mcthods appli ed. The gravit y fie lds 
of com plex 3D models, i.e. overhangi ng and non-overhangi ng sa[1 
domc~, are sim ulated. The gravimet ric anomalics for both cases have 
very simil ar shapes. Gravi ty modelling can d istingu ish between these, 
because the ex is ting mass di fferences result in anomal y d iffe rences 
bot h for surface profiles and X-sections. The capabil ity o f (he inver­
sion procedure is a lso show n in the d ise llssed synt hetic case. T he 
inversion manager is able to creale thc global struclUral forms repre­
sented "" a horiwn with constant density contrast (a two-l aycr model) 
from residual gravity anomal ies. 

L. INTRODUCTION 

Gravimetry is applied in hydrocarbon, groundwater, 
minera l and geotechnical studies. In hydrocarbon stud­
ies , grav it y modell ing is often used on a regional scale 
to deli ne3le [he basin and the main structural characte­
ristics assoc iated wi th the ta rgets. On a local scale it is 
used in struc tural inte rpre ta tion to de fine o il and gas 
struct ures. However, in de tailed reservo ir and produc­
tion modelli ng gravimetry is not fully ut ili sed , although 
a lot or d iffe rent da ta sc ts obtained by surface and bore­
ho le measurement s can be used to more precisely con­
strain Ihe interpreted structures (Fig. 1). 

Kljucne rijeci: gravi metrijsko modeliranje. kOlllpleksni 
modcl i, arhitektura softvera. 

Sazetak 
Razvi jena je clasli cna arhileklU ra so ft vera za gravi mctrijsko 

model iranje, Ie su izlozenc prednosli programa koji mogu obradivati 
kompt eks nc 3D mode lc. Soft vel' se sas toji od cetiri d ijela, koja su 
ugradcna U odgovaraj uce kompj utorsko succljc: programa za obliko ­
vanje i prikaz iva nje trodi mcnzio nal nih modc la (GOCA O soft ver, 
vcrzija 7.0), programa za fo rmaliranje i izmjenu podalaka izmed ll 
razlicitih so ftvera (GOCAD soflver iii GEOM O D .~oftver), programa 
s ugradenim al goritmima za racunanje g rav imc lrij skog djelovanja 
(Talwani -Ewing i Go tze- Lahmcycr metodc u pod ru cj ll metoda 
konacn ih clemenata) i prog rama Zll g rav imc trij sku inve rzij u tcme­
!jCllOg na Corde ll-Hendersonovoll1 postu pku inverzijc. 

Smat ra se bo ljom arhit ekturom drza ti odvojcn im so rlvc r za 
obl ikovanje i iz mj enu modcJ a od so ft ve ra za racunanjc g ravim elri­
jskog djelovanja, a posl upak inv erzij e sc os tvaruje razmjcnom 
podataka izmedu la dva dije!a. 

Na si nlc tskim pri mjerirna izloze na j e upo trcba i mog ucnos ti 
razv ij cnog son vera i primijcnjcnih metoda. Prika za na su teo re lska 
gravimel rijska polja solnih doma, kojc predslavljaj u siol.c llc gco loskc 
modele. Gravi melrijskc ano ma lijc ov ih modc la irnaj u vrl o s lic ne 
obrise, Ie je jedino gravimetrijski m modelir'lnjem moguce medu so­
bno razluCit i ove mod ele zbog razlic itog rasporcda masa, koje 
uzrokuju razliCi te vrijednosti anomalija, k:lko za pov['s inske profi le. 
tako i za X~sekcij u . Mog ucnos ti pos tu pka invcrzijc su prikazane na 
prilozenom si ntetskom primjeru. Pos tu pko rn inver.djc dobiven a je 
Slruk tu ra po gra nie i s razlikom gus toce, od nos no za jed nos tavni 
dvoslojni model, na temelj u rezidualnc gravimclrijske anomalijc. 

In order to use aU available data sets it is important 
10 have modellin g soft wa re whi eh sho uld be a t least 
three-dimcnsionallo accurately describe [he subsurface 
propcI1ies in space and enable realistic simulat ions to be 
form ul ated . It shou ld be ab le to process comp lex and 
deta iled three-dimensional models, as fo r surface mode­
ll ing. Therefore for borehole modelling, applying gravi­
metr ic a lgorithm s to boreho le measurements, can aid 
calculat ion of the e ffec ts o f both surface and boreho le 
local ions. The la tte r arc very impo rt an l in the produc­
tio n phase to constra in struc tura l inte rp rela tio ns of 
reservoi rs on a local scale. 

It must be emphasized that high-quality gravi ty data 
acqu ired by high-precision instruments shou ld be avail-

I University of Zagreb, Facul ty of Mining, Geology and Petroleum Engineering. Pi erollijev a 6, HR - l 0000 Zagreb. Croa! i:l. 

2 TNO Ins!itll ic of Applied Geoscience, Schoemakers traat 97, 2600 JA Del ft, Thc Nethcrlands. 



GRAVITY MODELLING , 
sasinmLD I 

Prospect models I 
• Reservoir modets I 

Production models . 

E XPLORATION P RODUCTION 

Fig. 1 Tile effec t of gravi t j modelling on buildill£ o f the <li llcrcll l 
geolog ical models ill prac tice from exploratiol110 product ion. 

able for M lC":CCSs rul g ravity modelling. On the local sca le 
hi gh-reso lut ion g ravi ty dat<l shou ld be used for g rav ity 
mode lling or s mail com plex st ructures . Data spac ing 
depends 011 the Illodel c le men t size, but it al so affects 
the precis ion or the reg ional -res idu<ll sepa ration, wh ic h 
is a vc ry c riti ca l s tep in gravi ty interpre tation fin d s ig­
nifican tly controls the final resu lts. 

Thi s pape r di sc ll sses a fl exible so ft ware archi tecture 
and the advantagc or having several programs to handl e 
complex three-d imensional modcls. 

2. G RA VlTY MODELLING 

T he proble m or g ravi ty interpreta ti on involves the 
solutions of bot h the forward and inverse problem. T he 
d irect problem or the forward modelling dea ls with the 
theore ti cal c al c u lat ion or the g ravi ty a no ma lies for a 
g iven model. The inverse problem o r g rav it y in vers ion 
deals wit h the calcu lation of structural fo rms fo r g iven 
g ravity a nomali es . Both prob lems can be solved by 
ll s ing two different approaches, th e s pace dornain 
approac h and th e freq uency domain o r spec tral 
appro,l(.;h. 

The "Space Domain Mcthods " ll sing a vo lulll e inte­
g ral , are frequently applied ror thc forward calcu lation 
of thc g ravitat iona l attraction of arbitrary shaped bod­
ies. In o the r te rm s, these arc ana lytica l techniques to 
cal c ul a te th e g ravity fields of complicated mode ls 
whose geome try is approximated by a certa in number 
or small e leme ntary bodies, such liS horizontal lam inas, 
prisms, or po lyhe drons . tha t have analytical so lu ti ons 
( Fig. 2). 

TA I. WAN I & EWING ( 1960) suggesled a method 
to cal culate the gravity anomal y of a three-dimensional 

a b 

body by the llullle rical integration or a sel of horizonta l 

polygona l laminae tha t approx imate the sh ape of the 
bod y. NAGY (1966) sta ted a formu la 10 ca lc ulate the 
gravitat ional allraction o f a ve rti cal rectangular p ri sm. 
All ite rative method for calculation of the g ravitationa l 
all rac ti on, whe re the causati ve body is approximated by 
ve rtical pri sms. was sugges led by CORDELL & IIEN­
DERSON (1968) . The gravity effcel o f each pri sm cle­
men t is approxim a ted by the vcni ca l right -cylinder 
sou rce rormu la and the ve rti cal - l ine-so urce ["orm ul a. 
BOTEZATU et a l. (1971) sugges led the 1l1elhod for cal­
culation of the gravi ty erfect by dividing the body into 
s1l1all-si7.ed cubes. Also M Uf'TI ( 1973, 1975) presenlcd 
approximation rormulas ror a cube's grav ity field and a 
g ravity modelling scheme whic h uses cubes as building 
blocks. The method based o n the exact ca lcula tions of 
the grav ity and magne ti c anomalies of po lygonal prisms 
was presented by P LOUr:r: ( 1976). BI-IATT AC fl A ­
RYY A (1978) considered computer mode lling in g ra vi­
ty and magnetic inte rpretat ion us ing the solu tion of the 
d irec t and inve rse problems. Three-dimensional in terac­
tive modelling in gravi ty anc! magne tics was sugges ted 
by GOTZE & LAHMEYER ( 1988). They presenLed the 
int e ractive modelling program where the body is con­
stru c ted from polyhedra of suitabl e geometry. The 
gravity c!l ec t of the homogeneous polyhedron is c~llcu ­

lated by tran sform ing the vo lume integral into a sum o r 
line integrals. 

The methods which lise Fourier nna lysis be long to a 
second g roup - " \Vavc Nu mbe r Domain Me thods". 
T hese methods conve rt the grav it y anomaly into Fouri ­
e r series and prov ide for g r<lvity in version, com monly 
us ing the Fas t Fourie r Trans form . Deve lo pme nt of 
these melhods sta rted w ith PARK ER ( 1973) and OLD­
ENBURGH (1974), who used Fourie r Trans fo rm Ana ­
lys is o f the gravity anoma ly. The re nre a large number 
of methods belonging to this g roup, but only a rew w ill 
be menLioned here. GRANSER ( 1987) useel Fast Fouri ­
e r Transform and an exponent ial dcnsi t y~dept h func­
tion. LINES ct al. ( 1988) ca rri e d out a coo perati ve 
invers ion 01 geophys ical dala. C II ENOT & DEBE­
GL IA ( 1990) presented the me th od or de pth -mapping 
invers ion or gravity ;:H1d magnetic fi e lds. GUSPI ( 1992) 
proposed one r rT gravity inve rs ion me thod based on 
Parkcr-Oldcn burgh formulas. 

Grav it y inversion us ing the rouri er transfo rm o r 
ot her techn iques leads to smoothing of the data because 
s uccessful convergence requires removal of the high­
frequency compone nts. Therefore it can be successfully 
<lp pl ied in the g ravi ty model ling o f large a re as (i.e. 
bas in s tructures) but not in the cases of small complex 

/ /:5\ 
/' \ 

/ ""-. 
c 

fi g. 2 Cross sec tiol1 or the booy 
represclHcd by vertical pri sms 
(a). cubic blocks (b), alld con­
tours (e). 



St1ll1anovac. RilSCm;l & Gril: System Archi tect ure for 3D Gravity Modelling 147 

INVERSION SCHEME MANAGER I 

network 

I 
MODEL BUILDERS I I DATA TRANSLATORS I I GRAVIMETRIC I 
AND VISUALIZERS AND TRANSMITTERS SIMULATORS 

structures (i.e. overhanging salt dome, overthrust, etc.). 
Gravity illversion can be very useful in generating an 
overview of the study area and thereby enabling separa­
tion o j' prospective zones and direction of further 
research activity. 

The forward modelling can be applied in the cases 
of small comp lex bodies and it allows more detailed 
interpretation of the gravity data. Besides, the inter­
preter can be more involved in the interpretation 
process applying a priori knowledge derived from th~ 
data acqu ired by other methods (seismic methods, the 
data from bore holes). 

3. ARCHITECTURE 

Different 3D models (velocity, geological, reser­
voir) that are produced by various model builders, arc 
not always accessible by gravimetric routines. Typical ­
ly, models have to be exchangcd and rcformatted before 
the algorithms will work. A good open architecture can 
keep the model building and updating software separate 
from the forward simulation software and can realize 
inversion schcmcs by communication between the two 
worlds. Thi s means that at any time, models of var ious 
re lated properties can be checked with different con­
strain ing data scts . It is also clear that the architecture 
applied is very Oexible enabli ng a relat ively rapid easy 
change of the implementations on the subsystems very 
quickly. 

The flexible architecture, consisting of four parts, 
can be implemented in a distributed compuler environ­
mcnt (Fig. 3): 

the three-dimensional model builders and visualizers, 

thc model representation translators, 

the forward simulation algorithms of grav im etric 
data, 

the inversion (model updating) scheme manager. 

3.1. 3D MODEL BUILDERS AND VISUALIZERS 

Two representations are typically used in 3D model 
builder programs, i.e. grid or finite clement representa­
tions (FEM) and boundary representations (BREP), 
each of them has advantages ror specific visualization 
and simulation algorithms. In this work, the connect ion 
has been established bet ween a three-dimensional 

Fig. 3 Software architecture . 

BREP modcllcr (GOCAD) and several gravimetric pro~ 
grams operating on different FEM representation s of 
the three-dimensional model. 

GOCAD is three-dimensional modelling software 
based on the "Discretc Smooth Interpolation Method" 
which is the centre of a new approach based on decoIll­
position of lines into segments, surface into triangular 
facets and volume into tetrahedrons (MALLET, 1989). 
GOCAD can be used very efficiently for modelling sur­
faces dcl'ined by various type data, and for exporting 
the data for the calculation processes. In the export 
process a transformation is needed in some cases. It is 
especially important that GOCAD can eas ily handl e 
complex structures such as overhanging salt domes and 
overthrusts. It also enables the easy calculation of syn­
thetic gravimetric data for such geological models. 

3.2. MODEL REPRESENTATION TRANSLATORS 

Several transformation functions were needed and 
implemented in the dist ributed computer environment 
to exchange data from the model builder to the different 
gravimetric simu lation routines. The internal GOCAD 
modeJ representation, i.e. triangulated surface (BREP) 
was transformed into two finite clement (FEM) repre­
sentations, being horizontal contouriine slabs, and verti­
cal triangular prisms. This is realized using avai lable 
transformation programs in GOCAD or external pro­
grams such as GEOMOD. The output is used by for­
ward and inverse modelling software. 

3.3. FORWARD GRAV IMETRIC SIMULATION 
ALGORITHMS 

Two forward calculation algorithms have been 
implemented, Talwani-Ewing (using horizontal slabs) 
and Gotze-Lahmeycr (using vert ical prisms). 

In the Talwani-Ewing method the three-d imensional 
body is first represented by contours that arc then 
replaced by horizontal irregular n-sided polygonal lami­
nae. If the number of sides n is suffic ien tly large , the 
polygon can approximate the contour line as closely as 
desired. The gravi ty anomaly produced by each lamina 
can be calculated analytically at any external point. A 
contour line on the surface of the arbitrary shaped body 
at depth z below stat ion P (Fig. 4) is replaced by an n­
sided polygonal lamin a of infinitesimal thickness dz 
(T ALW ANI & EWING, 1960). The gravity anomaly at 



148 

y 

x 
p ,L---------+ 

-------
z 

y 
H 

G 

F 

contour at depth z 

Fig.4 Geometri c elements of the polygonal lamina used 10 calculate 
the gravi ty effect (a fter TALWANI & EWING, [960). 

station p. whic h is the ve rti cal compo nent of the g rav i­
tat ional a ttrac tion , prod uced by me ntioned polygonal 
lamina can be expressed as: 

g, = V dz 

( I ) 

Te rm V is Ihe g ravi ty anomaly caused by the poly­
gona l lamina pe r unit thi ckness. V is expressed by the 
surface integral, but it can be reduced on two line inte­

grals (both along the boundary of the polygonal lam ina) 
and can be g iven by the expression: 

(2) 

where K is the un iversa l gravitational cons tant , p is the 

volu me de ns it y o r the la mina a nd z , y a nd r arc the 
cy l indrical coordinates (Fig. 4). P ' is the projection of 
thc station P on the plane of the polygonal lam ina. 

A ft e r severa l transformations and substi tutions the 
fin a l expression for a gravity anomaly caused by poly­
gonal lam ina is: 

zcosil , }] 

) Pi2 + Z 2 

(3) 

GeOlogia Croalica 49(2 

x 

y 

Fig. 5 Geomelry o f the etementary body used to calcu lme the gravity 
CrrCCl (after GOTZE & LAHME YER, 1988). 

whe re Pi:::::PllJ. All terms (Pi' Yi' Yi+i' cosch and cosy) can 
be e xpressed in Xi and Yi coordi na tes of the polygonal 
points. 

The tota l anomaly g, produced by who le body is 

expressed as: 

'-
g, tota l = f dz 

'-
(4) 

Us ing the afore me ntioned express ions the aut hors 
gave the ex pression sui tab le fo r programmi ng. w he re 
a ll terms arc given in x , y a nd z coo rdina tes of the 
polygonal points. 

The Gotze-Lahmeyer method to calculate the allrac­
tion of a homogeneous polyhedron al station P is based 
on the eval ua tion of the vo lume integra l. The gravity 
anomaly g(P), after derivation of the pot e nlia l with 
respect 10 the vertica l component z, is g iven by exprcs­
s ion: 

au fff a 1 Tz (P) = g(P) = Kp az (e) dv 
Poly. 

(5) 

w here U(P) is the potential al s tat ion P, r is d is tance 
be tween P a nd dm (dm=pdv, clv= pdxclydz), a nd K is 
the gra vi tationa l cons tan t. Correspo nd ing to the rele­
valli theories of vec tor ana lys is, a n equation is estab­
li shed where the surface integral has to be calculat ed 
fo r the w hole polyhedron sur face. A transformati on o f 
the coordinates syste m (into the surface oriented coor­
dinate system) is carried out to s implify the mathemmi­
cal expressions in the evalual ion of the surface integral. 
Finally a conversion of the surface integra l into a linear 
in tegra l via polygon, which limits the surface, is carried 
out. 



SUlll<l I IUV;I(:. RIISC111<l & I3ril: SystCl ll Architecture for 3D GmvilY MlXlclling 14 9 

Fig. 6 Overhanging salt dome Wilh gravity field at surface (a), which is detached as a map (b) . Contour interval 0 .2 x I05m/s2, density CO!l1rast 
t.p=-400 kg/m). 

The au thors also produced a final formula to calcLl­
late the grav ity anomaly produced by the polyhedron , 
wh ich is suitable for programming. A subroutine has 
been published, written in FORTRAN, 10 calculate Ihe 
gravity effect of a vertical prism , infinite in the z direc­
tion with a triangle as the upper boundary. This elemen­
tary body, whose shape (geometry) is presented in Fig. 
5, ean be used to construct more complicated bodies. In 
th is case the gravity effect of the entire body would be 
expressed as a sum of the gravity effects of the elemen­
tary bodies. This s ubroutine was the second routine 
implemented for the gravity forward modell ing. 

The option, which uses the Gbtze-Lahmeyer method 
can be used to calculate the gravity effect of noncom­
plex (sim ple) structures. For computation of the gravity 
effect or complex structures, as for instance a salt dome 
(overhanging) or overthrust, the option using the Tal­
wan i-Ewing method should be used. 

3.4. INVERSION SCHEME MANAGER 

The inverse sc heme manager should support the 
inversion strategies. One approach is by a trial and error 
process, another by automatic (constrained) inversion. 

The Cordell-Henderson gravity inversion scheme is 
use rul to producc the first model, which is constrained 
here on a simp le two-layer modeL It is an iterative 
scheme based on the standard trial and error approach. 
The scheme consists in general, of threc aspects: 

construction of a simple start model , 

calculation of the gravity anomaly for a given model, 
and 

modifi cation of the model based on differences 
between calculated gravity anomaly and observed 
gravi ty anomaly. 

In the first step the gravity anomaly is digitized on a 
rec tangular grid . The interpreted body is represcnted by 
rectangular vertical prisms in a regular grid with a f1aL 
top and bottom. The thickness of the prism elcment 
below the q-th grid point (station) tl.q for the starter 
model is es tabli shed by means of the Bougucr slab for­
mula: 

t 1.q = BgobS,q 8=_1_ 
2nKp 

(6) 



ISO Geolog ia Cro;llica 49n 

Fig. 7 Cylindrical ~alt dome with gravity field ,II surface (a), which is dClached as a 1ll;IP (b). Contour inlcrval 0.2 x I0 5m/s2, dcm;il y contrast 
6p=.400 kg/mJ

. 

where K is the g ravitat ional constan t, p is the density, 
goo.,q is the observed gravi ty at q-th station. The authors 
sugges ted the following relationsh ip for the next modi ­
fi cation or improvement of the model: 

t t (g".,q) 
n+I ,q = n,q g-­

calc ,n,q 

(7) 

where J1 is number of iterat ion and gc~lc.n.4 is calculated 
gravity value below q-Ih station at n-th iteration. 

The prisms are thcn approximated by vertical cyl in­
ders si nce the fo rmu la calcu lating the grav ity e ffect o f 
the ven ical cylindcr is much simpler and consumcs less 
calculati on tim e. The goodness of fit is expressed in 
terms of the root-mean -square error (RMS error). How­
ever, also the largest error can be used. 

The inversion scherne described is used for obtain­
ing the initial model of a two-layer structure. The result 
or more complex models (based on ot her da ta) can be 
verifi ed by using fo rward modelli ng techniques <.I nd 
adjusted by Ll sing the trial and error approach. 

4. CASES 

Two synthetic cases (o verhanging and non-over­
hanging salt dome) are shown as an example or thc use 
or architecture and the capability to simulate the gravi­
mct ric fields fo r complex mode ls (figs. 6 and 7). The 
mode ls and the gravimet ric fields a t the surface arc 
shown us ing GOCAD software ror bO lh cases. The 
erfect of the difference in shape of the bodies is clearly 
visible (Figs. 6 and 7). The g ravity anomalies at the 
cent res of the sal t domes arc -3 .6 x 10-5 111 / S2 for the 
overhangi ng and -4.6 x 10·5m/52 fo r the cy lindrical 
dome, respect ively. 

Both surface gravimctric anomalies have conccllI ric 
shapcs, but on the surface gravity profi les (Fig. 8) the 
differences in gravity attraction arc ev ide nt , a lthough 
they decrease with the increase or the distance [rom the 
centre or domes. Calculated gravimet ric data 1'01' X-sec­
tions just outside the cy lindrica l and overhanging salt 
domes arc shown in Figs. 9 and 10. Fo r cylindrical 
body, the isoli nes arc almost ve rl ica l for th e fi rst 500 
metres, and fo r overhang in g structure, thc gradient 
increases faste r wi th dcpth. The X-section of d iffe ­
rences in grav ity attract ion bet ween dom cs arc shown 



.sumanov;lC. Ritsema & I3ril: System ArchilCCllIrc for 3D Gravity Modelling 

o 

-1 

'iii 
"' g 

-2 
~ 

E 
0 
c 
0 

-3 
~ 

:5 
~ 
rn 

-4 

~p=~400 kglm3 

-5 
-1500 ~1000 ~500 0 500 1000 1500 

distance from centre of dome (m) 

Fig. ~ Gravity anolllalies al the surface profile. 

in Fig. 11, where the values also decrease with the 
increase of distance from the centre of domes. 

With respect to the previous examples it can be con­
cluded that it is impossible on the basis of anomaly 
shape to make any conclusions about source body. But 
gravity modelling can distinguish between cases 
because the existing mass differences result in clear 
anomaly dirferences both for surface profiles and for X­
sections, although the possibility of the prccisc separa­
tion of the domes decreases with the increasing distance 
from thc dome ccntrc. Thcrefore in a practical sense, a 
salt dome can be defined more precisely if bore hole 
and gravity measurements afC situated closer to the 
domc ccntre. 

The best illustration of the previous statemel1l is the 
following example. If automatic gravity inversion is 
applied on the gravity field produced by an overhang­
ing dome, satisfactory results can not be obtained. The 
automatic gravity inversion procedure can not process 
complex structure, which has few z-values for the same 

surface 

·200 

·400 

I 
~ ·600 

\ :-li ·0 • 
" -800 

·JOO!) ~"~ 
·1200 

0 2(1) 41)0 600 8110 10(1) 12110 1400 

distance from centre of dome (m) 

Fig. 10 Gravity field for X-sect ion of overhanging salt dome. Con­
lOur interval 0.5 x 1O-5m/s'. density contrast ~p=-300 kg/m3. 

lSI 

surface 

-2flfl 

-400 

I 
~ ·61JO -
li 
m 

" 
·800 

~o " .6' 

·1000 

·1200 L....--:.--"--c
2
"00::--"C

400
:.,-

600 800 1000 

distance from centre of dome (m) 

rig. <) Gravity field for X-sect ion of cylindricdl salt dome. Contour 
interval 0.5 x 1O-5m/s2

, density contrast ~p=-300 kg/mJ . 

x and y coordinatcs. The result is unstable inversion, 
which mcans that a satisfactory model cannot be pro­
duced. The conclusion is that such a structure has a 
complex shape, which can be established by forward 
modelling. It is obvious, that gravity modelling can dis­
tinguish these cases. Other complex structures, sLlch as 
rcservoirs fissurcd by faults or other structural features 
will also benefit from these kinds of analysis . 

Anothcr synthetic case can be described to illustrate 
the capabi lity of the inversion procedure. The inversion 
manager gives the structural form s of the horizon with 
constant density contrast bounding a single laycr from 
residual gravi ty anomalies using a regular grid . It must 
be noticed that the res id ual-regional separation should 
be carried ou t, mcaning that the regional gravity effect 
must be removcd from obscrvcd gravity anomalies by 
using appropriate methods. The gravity ficld at the sur­
facc produces the structure shown in Fig. 12. Satisfac­
tory interprctation of the gravity field can be attained in 
thc casc or morc comp lex geologica l models , as well. 

surface 
0 

.100 

.4(1) 

I 
'K ·600 

m 

" 
-800 

·1000 

·1200 
200 400 1200 1400 

distance from centre of dome (m) 

Fi g. 11 Di fference in gravity fields of overhangi ng and cylindrical 
sal t domes. Contour interval 0.1 x 1O-5m/s' . 



152 Geologia Croatica 49(2 

Fig. 12 The slructure obtained from a given gravily field (gravity inversion). 

5. CONCLUSIONS 

Two conclusions can be drawn. Firstly, the applica­
tion of gravimetric data to constrain 3D simple and 
complex (structural and reservoir) models in the appra­
isal and development phases is possible and should be 
incorporated in daily practice. Secondly, the described 
architecture, separating the model building and visuali­
zation subsystems from the algorithmic subsystems for 
forward and inverse simulation will support this in an 
open environment. The communication between two 
environments consists of a set o[ transformation rou ­
tines for various three-dimensional model representa­
tions. It will allow the easy intcgration of other simula­
tion software, such as magnetic, seismic and so 011. 

Acknowledgement 

We thank Prof. Dr. Zeljko ZAGORAC for reading 
of the manuscript and useful comments. We also wish 
to thank Prof. Dr. C. REEVES (ITC - International 
Institute [or Aerospace Survey and Earth Sciences, 
Delft, The Netherlands), who acceptcd the manuscript 
for review, and Dr. Sally BARRITT (lTC, Delft) for 
critical review of the manuscript and suggestions. 

6, REFERENCES 

BHATTACHARYYA, B.K. (1978): Computer mode­
ling in gravity and magnetic intcrprctation. - Geo­
physics. 43, 912-929. 

BOTEZATU, R., VISARION, M., SCURTU, F. & 
CUCU, G. (1971): Approximation of the gravita­
tional attraction of gcological bodies. - Geophysical 
Prospecting, 19,218-227. 

CHENOT, D. & DEBEGLIA, N. (1990): Three-dimen­
sional gravity or magnetic constraincd depth inver­
sion with lateral and vertical variation of contrasl.­
Geophysics, 55 , 327-335. 

CORDELL, L. & HENDERSON, R.B. (1968): Iterative 
three-dimensional solution of gravity anomaly data 
using a digital cOIllputer.- Geophysics, 33, 596-60l. 

GOTZE, Il. -J. & LAHMEYER, B. (1988): Application 
of three-dimensional intcractive modeJing in gravity 
and magnetics.- Geophysics, 53, 1096-1108. 

GRANSER, H. (1987): Three-dimensional inlerpreta­
tion of gravity data from scdimcntary basins using 
an exponential density-depth fUllclion.- Geophysica! 
Prospecting, 35 , 1030-1041. 

GUSPI (1992): Three-dimensional Fourier gravity 
inversion with arbitrary density contrast.- Gco­
physics, 57 , 131-135. 



SlIl1l;mOV<ll:, Ril ~elll:l. & I3ril: System Architecture ror 3D Gravity Ivl odcll illg 

LINES , L.R. , SCHULTZ, A.K. & TREITEL, S. 
(1988): Coopcrative invers ion of geophysical data.­
Gcophysics, 53, 8-20. 

MALLET, J. L. ( 1989): Discrete smooth intcrpolation. ­
ACM-Transactions on Graphics, 8/2, 121 - 144 . 

MUFTI , I.R. (1973): Rapid determination of cube 's 
gravity ficld.- Ge ophy sical Prospecting, 21, 724-
735 . 

MUFTI , 1.R. (1975): It erative gravity mode lin g by 
us ing cubical blocks.- Geophysical Prospecting, 23, 
163-198. 

NAGY , D. ( 1966): The grav itational attraction of a 
ri ght rectangular prism. - Geophysics , 31, 362-37 1. 

153 

OLDENBURGH, D.W. (1974): The invcrsion and 
intcrpretation of gravity anornali cs .- Geophysics, 
39, 526-536. 

PARKER, R.L. (1973): Thc rapid calculation of poten­
tial anomalies. - Geophy s. J. Roy . A s t!". Soc. , 41, 
447-455 . 

PLO UFF, D. (1976): Gravity and magnctic fi clds of 
polygonal prisms and application to magnetic ter­
rain corrections.- Geophysics , 41 ,727-74 1. 

TALWANl , M. & EWING, M. (1960): Rapid computa­
tion of gravitational attract ion of thrce-dimensional 
bodies of arbitrary shape. - Geophysics, 25, 203 -225. 

Manuscript received January 11 , 1996. 
Reviscd manuscript accepted September 11 , 1996. 



Geologii1 Croatic~ 49/2