1 Introduction

Several fascinations of the ancient Egyptian engineering exist today, one of these fascinations are obelisks. A good number of obelisks are currentlylocated in major cities in Europe and North America after being transported fromEgypt, in addition to other obelisks still existing in Egypt. The number of obeliskswith a height exceeding 10 m could be more than 50 [1]. Natural disasters and soil problems areconsidered the main reasons for the drop in numbers of obelisks existing todaycompared to thousands of years ago. One of the most famous obelisks is the onecurrently located in the Vatican. This obelisk was moved to Italy during the Romanrule of Egypt and then repositioned to the Piazza di San Pietro in the sixteenthcentury. On the other hand, the tallest of all ancient Egyptian obelisks is theLateran obelisk. That obelisk was raised in Laterano in Rome in the sixteenthcentury [1].

Meanwhile, historical and geological evidences confirmed that ancient Egyptian obelisks were carved from the granite formation located in Luxor, Egypt[2]. The granite from which theobelisks are carved is generally called the “Red Aswan Granite”consisting of reddish feldspar crystals together with quartz, plagioclase andbiotite [3] [4]. The properties of this granite have beenexperimented by [5] and it was foundthat this material had an average compressive strength of 140 MPa which isnearly quadruple that of concrete however it experienced a brittle mode of failurewith no plastic deformation experienced by any of the tested samples. Within thatsame study, it was found that the average modulus of elasticity was 5.4 GPa.[6] studied the response of fivedifferent obelisks when subject to seismic loading. That was done by performing afree vibration analysis followed by a time-history dynamic analysis. The highmodulus of elasticity of this material (when compared to other natural stones) wasfound to be the main reason for the relatively low lateral deflections and lowstresses when subject to earthquakes [6].

Furthermore, one of the main things that makes obelisks unique from awind engineering perspective is their slenderness. The height to width ratio of suchstructures could vary from 9 to 12 which makes such structures considered to besignificantly slender. This is a point of similarity that these structures sharewith tall buildings as a sky-scraper could approach that range of slenderness ratiohowever the difference in scale between sky-scrapers and obelisks is expected tocause a different response to wind loads. Extensive research was conducted on theresponse of tall buildings to wind loads by numerous researchers like [7] [8], [9], [10] and [11] within the past four decades however nobody has studied thebehavior of obelisks under wind loading.

Most of the research studying the pressure coefficients of tallbuildings was either using field results that are very difficult to acquire orrelatively expensive wind tunnel testing or CFD modeling. It has been proven thatalthough field results and experimental results are the most accurate techniques tomodel tall buildings subject to boundary layer winds, CFD modeling could providesufficiently accurate results especially with proper meshing and boundary conditions[12] [13] [14]. However, the slenderness of the buildings studied by thevast majority of researchers was less than the slenderness of a typical ancientEgyptian obelisk which could have a height-to-width ratio ranging from 9 to 12. Thisfact could cause the wind pressure distribution on the surface of an obelisk to besignificantly different than that of a tall building in addition to the fact that itwill make it less stiff.

Another unique feature of ancient Egyptian obelisks is the taper thatthey have. The cross section dimension of each obelisk varies with its height with asignificant taper that could have a taper ratio ranging from 26 to 43 [1]. This taper may affect the pressuredistribution on the surface of each obelisk. Some of the research studying theresponse of tall buildings to boundary layer wind was performed on buildings with ataper however not as that of ancient Egyptian obelisks and not as slender[15]. Additionally, no researcherso far has studied the effect of varying the taper of any structure on the pressurecoefficients.

Meanwhile, another interesting feature about obelisks is that althoughthey do have a taper, they always maintain a squared cross-section causing each twosuccessive translational modes of vibration to have identical natural periods due tohaving the same stiffness and mass distribution within the two perpendicularhorizontal axes [6]. This is also thereason why the torsional mode of vibration was a very high mode causing no torsionalcomponent of the response when such structures were subject to lateral earthquakeloads as the perfect symmetry caused the center of mass and the center of rigidityto coincide negating the presence of any torsional effects [6]. A similar behavior is expected to happen ifwind loads are applied on such structures provided that these loads are symmetric onthe structure itself; however, no researcher has studied the response of ancientEgyptian obelisks to wind loads till now further than studying the effect of taperand height variation.

Unfortunately, until this research at hand has been initiated, noresearch was done to study the variation of the wind pressure and the pressurecoefficient on the surfaces of slender structures with a taper similar to that ofancient Egyptian obelisks. Furthermore, nobody studied the structural response ofthese obelisks under any type of wind loading to determine whether such structureswill fail due to wind loads or not. This scarcity of information regarding theresponse of obelisks to wind loading initiated the need for the current presentedresearch.

The objective of this research is to study the pressure variationwithin obelisks when subject to boundary layer wind with the obelisks’heights and taper and under different angles of attack (θ). Furthermore, itis necessary to determine the most critical wind load case that will cause thelargest stresses within each of these structures. This is necessary to determinewhether such structures could fail due to wind loads or not.

In order to achieve this goal, CFD models were prepared for two of thelongest existing ancient Egyptian obelisks located in Lateran and the Vatican havingdifferent heights, cross-sections and tapers as shown in Fig. 1. The CFD modeling was performed for each obeliskunder four different angles of projection to the wind load as shown inFig. 2 and for three differentreference velocities hence performing twenty four different CFD models. The outputof this process was used to perform structural analysis for each of the studiedobelisks under their different load cases representing the different angles ofattack. Finally, the most critical load case was determined for each obelisk basedon the maximum stresses identified from the structural analysis stage. The endresult is a quantification of the variation in the pressure coefficient on thesurfaces of the two ancient Egyptian obelisks under study. Additionally, thestresses within the two obelisks under different angles of projection to boundarylayer wind are acquired, providing an assessment of the obelisks’ structuralperformance that will inform the scientific community whether these obelisks need tobe further protected from wind loads or not.

Fig. 1
figure 1

(a) Lateran obelisk (b) Vatican obelisk

Front and top views of Lateran and Vaticanobelisks

Full size image
Fig. 2
figure 2

The four different studied loading directions

Full size image

2 Materials and methods

2.1 Mathematical modeling

The present mathematical research was built on solving 3-Dprincipal equations that described airflow around two ancient Egyptian obelisks(one at a time) by the CFD program fluent ver. 19.0 [17]. This mathematical approach solves thepartial differential equations (PDEs) governing the transport of mass, threemomentum, in a fully turbulent three dimensional domain under steady state andincompressible ideal gas conditions. In addition the standard k–εmodel equations for turbulence closure were used. The different governingpartial differential equations are usually expressed in a general form as:

$$\frac{\partial }{\partial x}\rho U\phi +\frac{\partial }{\partial y}\rho V\phi +\frac{\partial }{\partial z}\rho W\phi =\frac{\partial }{\partial x}\left({\varGamma }_{\phi }\frac{\partial \phi }{\partial x}\right)+\frac{\partial }{\partial y}\left({\varGamma }_{\phi }\frac{\partial \phi }{\partial y}\right)+\frac{\partial }{\partial z}\left({\varGamma }_{\phi }\frac{\partial \phi }{\partial z}\right)+{S}_{\phi }$$
(1)

Where ρ represents the air density, \(\phi\) represents the dependent variable, \({S}_{\phi }\) represents the source term of\(\phi\), (U, V, W) are the velocity vectors, and\({\varGamma }_{\phi }\) represents the effective diffusioncoefficient.

The standard k-ε model was chosen based on an earliercomprehensive verification study of the buildings for theaerodynamics of a building, including the standard, RenormalizationGroup (RNG) k-ε model and realizable, the standard k-ωmodel, Large Eddy Simulation and the Shear Stress Transport (SST). Astudy, presented by [16], showed that the k-ε model was capable ofmore accurately predicting aerodynamic clouds, with a variation thatis less than 3% compared to the result of the correspondingwind tunnel.

2.1.1 Mesh generation

Figure 3represents an ancient Egyptian obelisk, which has main domain dimensions(12 H × 3 H × 4 H) (length ×height × width) where H is the obelisk height, meshed with more than1.936 × 107 tetrahedralcells. The mesh on the obelisk wall has a size of 0.1 m with aninflation on the obelisk walls. The first layer value is 0.1 m andincreases to 1 m with the remaining whole domain, and the rest of thefield is 1 m. The domain size has been optimised, to optimize themesh size and the computing time. This control volume has been previouslyproven to be sufficient [18].Fine meshing around obelisk walls is used to precisely capture boundarylayer characteristics and hence boost model dependability. A gridindependence research was carried out to guarantee the numericalsolution’s stability and convergence, as well as its independencefrom the mesh size chosen. The exterior model of the configuration iscreated using well-known solid modeller Pro-Engineer and ANSYS ICEM-CFD. TheICEM-CFD generated surface configuration is used to generate a good mesh inthe form of a mesh file. Fluent ver. 19.0 solver is used to solve thegenerated domain matrix equations while CFD-post is used to imagine theresults. Figure 4demonstrates grid independency check. The mesh dependency was studied bysolving the flow field for seven mesh configurations made of(2.35 × 107,1.936 × 107,1.6742 × 107,1.3678 × 107,1.0968 × 107,8.825 × 106 and6.621 × 106 cells)respectively, results presented that 8% variance in the pressurecoefficient between the coarser and finer mesh and less than 0.07%variance exists between the two finer meshes.

Fig. 3
figure 3

Actual geometrical configurationinvestigated

Full size image
Fig. 4
figure 4

Maximum pressure coefficient at a diverse number ofmesh cells(θ = 45o)and 20 m/s for Lateran obelisk

Full size image

2.1.2 Boundary conditions

At the inlet, uniform velocities of 20, 25 and 30 m/s(with different angles of projection to wind 0o,15o, 30o and45o) are applied with a turbulence intensityof 5%, demonstrating air movement due to wind speed. When based onthe building heights, the corresponding Reynolds numbersRe = ρUH∕µ for Lateran and Vaticanobelisks had values of(40.24 × 106 –47.44 × 106) and(31.62 × 106− 40.24 × 106),respectively. Based on the building width, the Reynolds numbers are(3.74 × 106 –5.6 × 106) for Lateranand (3.35 × 106− 5 × 106)for Vatican obelisks [19]. Theobelisk body surface was demonstrated as a no-slip boundary wall with zeroroughness representing the actual condition of the obelisks due to thenature of the granite surface from which they were made. For the lowest andupper side’s boundaries of the domain, a slip-wall boundary(symmetry) was used. At the exit of the computational domain, atmosphericstatic pressure was selected.

2.1.3 Validation

The CFD model used within the current study was compared to apreviously published research paper that provided an experimental andnumerical comparison of mean pressure coefficient Cp within the wind tunnel[13]. This comparison wasperformed to examine the validity and applicability of the current CFDmodel. The standard k-ε turbulence model is employed in thesimulation and the results are compared against the experimental results[13], as presented inFigure ‎5. This model waspicked from a pool of several different turbulence models because it matchesthe experimental data well.

Fig. 5
figure 5

Standard k-ε Simulation Results Compared toWind Tunnel Measurements [13]

Full size image

2.1.4 Numerical Solutions

Pressure-velocity coupling was selected as the coupledalgorithm, pressure interpolation was PRESTO and 2nd -order upwind finitedifference schemes were used for both the viscous terms and the convectionterms of the governing equations. Gradients are computed with thegreen-gauss node based method. The simulations were achieved with thecommercial CFD code ANSYS Fluent 19.0, which uses the control volume method.Convergence was carefully checked and the iterations were finished when notthe entire residue showed any additional decrease with an increase in thenumber of iterations. Besides, the scaled residuals were about10− 4 for continuity, turbulentkinetic energy, turbulence dissipation rate and10− 5 for momentum. Whilesolving the program, important variables like (Cp)should be monitored to ensure that the solution is convergent.

The pressure coefficient is defined as:

$${C}_{p}=\frac{(P-{P}_{\infty })}{0.5 \rho {{U}_{\infty }}^{2}}$$
(2)

Where the reference static pressure (atmospheric pressure)(Pa), P is the static pressure(Pa), ρ is the density of air(kg/m3) and is uniform inlet velocity(m/s).

2.2 Structural analysis

The four different angles of attacks studied represent the fourdifferent wind load cases studied for each of the two obelisks. The pressure wasmultiplied by the surface area to give the load perpendicular to the surface.This load was applied within the structural analysis process to acquire thestresses within each obelisk under each of the four studied load cases.

Each of the two obelisks was modeled on SAP2000 [20], in which 3-D 8-node solid finiteelements were utilized. The choice of such elements was to represent thestiffness along the obelisk height in the most accurate way as the significanttaper within each obelisk will cause a variation within its moment of inertiaalong its height. If 1-D or 2-D elements were used, the analysis would have beensignificantly in-accurate within at least one dimension; however, using 3-Delements guarantees the highest achievable accuracy in the analysis. Each of thetwo obelisks was modeled in its original dimensions representing the conditionof each during the era in which it was originally carved.

Each of the four load cases representing the different angles ofattack has been applied on each of the two modeled obelisks. The forces producedby the CFD were used in the structural analysis of each of the two obelisks whensubject to a 30 m/s wind event at the four different projection angles.The target of this exercise is to determine the most critical wind load case foreach of the two obelisks.

In order to determine whether a static analysis is sufficient or adynamic analysis will be needed, the natural frequencies of the obelisks werecompared to the dominant frequencies of winds. The natural frequencies ofvibration of the two modeled obelisks are 1.015 Hz for the LateranObelisk and 1.308 Hz for the Vatican Obelisks as reported by[6]. These values of naturalfrequencies are significantly higher than the values of dominant frequencies ofwinds which could be as low as10− 3 Hz or even10− 4 Hz [21] [22] [23]. Hence,the analysis performed was a static analysis as there was no need to analyzethese structures dynamically.

3 Results and discussion

The pressure coefficient contours at 60% of the height of theVatican obelisk are shown in Fig. 6.As expected, a symmetric distribution of pressure is shown inFig. 6a and d representing the0o and 45o casesrespectively as each of these two cases is a symmetric load case and thedistribution of pressure is expected to be symmetric for these cases. Meanwhile, thetwo asymmetric load cases of the 15o and30o have caused an expected asymmetric pressuredistribution as shown in Fig. 6b andc. A very similar behaviour is seen when examining the pressure coefficient contoursat 60% of the height of the Lateran obelisk are shown inFig. 7 in which a symmetricdistribution of pressure for the symmetric load cases is shown inFig. 7a and d while, the twoasymmetric load cases have caused an asymmetric pressure distribution as shown inFig. 7b and c. However, althoughthe distributions of pressures for the two obelisks look similar, the values of thepressure coefficients for the Vatican obelisk are different from their counterpartsin case of the Lateran obelisk. For most of the cases, the values of the pressurecoefficients for the Vatican obelisk were higher than their counterparts in theLateran obelisk. This could be attributed to the differences in taper andslenderness between the two obelisks as the Lateran obelisk is more slender than theVatican obelisk while the Vatican obelisk has a more significant taper as shown inFig. 1.

Fig. 6
figure 6

The pressure coefficient contours at 60% height ofthe Vatican Obelisk for 4 different angles of windprojection

Full size image
Fig. 7
figure 7

The pressure coefficient contours at 60% height ofthe Lateran Obelisk for 4 different angles of windprojection

Full size image

Furthermore, the variation of the maximum positive pressure coefficient(Cpmax) along the height of the Vatican obelisk andLateran obelisk is shown in Fig. 8.The maximum positive pressure coefficients for the Vatican obelisk shown inFig. 8 are slightly larger thanthe maximum positive pressure coefficients for the Lateran obelisk. However, itcould be noticed that for each projection angle the differences inCpmax between the two obelisks are negligible. B thegeneral trend of a nonlinear increase of the maximum positive pressure coefficientwith the height very similar to the boundary wind velocity profile is common betweenthe two obelisks as shown in Fig. 8.

Fig. 8
figure 8

Variation of the maximum positive pressure coefficient alongthe height of the two obelisks for 4 different angles ofprojection

Full size image
Fig. 9
figure 9

Variation of the maximum negative pressure coefficient alongthe height of the two obelisks for 4 different angles ofprojection

Full size image

Meanwhile, the situation seems different for the maximum negativepressure coefficient (Cpmin) that is shown inFig. 9 as a unique trend wasexperienced for each of the four angles of wind projection studied. In general, andfor both obelisks, the cases with small angles of 0o and15o have experienced higher magnitudes of suctionthan the larger angles of 30o and45o. However, the trends of the variation of themaximum negative pressure coefficients of the Vatican obelisk shown inFig. 9 are different from thoseof the maximum negative pressure coefficients of the Lateran obelisk. Thesedifferences are more significant when examining the cases with small angles of0o and 15o. Thissignificant change could be attributed to the fact that the suction for these caseshad higher values than the other two cases with larger angles and these highervalues are more sensitive to the differences in slenderness and taper than the lowervalues experienced by the 30o and45o angles as the taper and slenderness of these twoobelisks are significantly different as shown in Fig. 1.

Figures 10 and11 show the vertical variation of thepressure coefficient contours at an angle of attack of 0ofor the Vatican and Lateran obelisks respectively. The pressure coefficient contoursin the vertical plane perpendicular to the wind are shown inFigs. 10a and 11a while those in the vertical plane parallel tothe wind are shown in Figs. 10b and11b. The variation in the pressurecoefficients in the vertical plane perpendicular to the wind are very similar forboth obelisks showing an increase in the negative pressure with height till reachingthe maximum negative pressure at the location of the pyramidon at the top of eachobelisk however the difference in taper and difference in slenderness between thetwo obelisks caused the negative pressure on the Vatican obelisk to be higher thanthat acting on the Lateran obelisk as the Lateran obelisk has a higher taper and ismore slender as shown in Figs. 10aand 11a.

Fig. 10
figure 10

Vertical variation of the pressure coefficient of theVatican Obelisk for a 0o angle of windprojection

Full size image
Fig. 11
figure 11

Vertical variation of the pressure coefficient of theLateran Obelisk for a 0o angle of windprojection

Full size image

On the other hand, the variation in the pressure coefficients in thevertical plane parallel to the wind are very similar for both obelisks showing anincrease in the positive pressure with height in the windward direction. However,the difference in the taper of the pyramidon in each obelisk caused the pressurecoefficients to significantly differ between the two obelisks at the location of thepyramidon as shown in Figs. 10b and11b.

Fig. 12
figure 12

Vertical variation of the pressure coefficient of theVatican Obelisk for a 45o angle of windprojection

Full size image
Fig. 13
figure 13

Vertical variation of the pressure coefficient of theLateran Obelisk for a 45o angle of windprojection

Full size image

Figures 12 and13 show the vertical variation of thepressure coefficient contours at an angle of attack of45o for the Vatican and Lateran obelisksrespectively. The pressure coefficient contours in the vertical plane perpendicularto the wind are shown in Figs. 12aand 13a while those in the vertical planeparallel to the wind are shown in Figs. 12b and 13b. Thevariation in the pressure coefficients in the vertical plane perpendicular to thewind are very similar for both obelisks showing an increase in the negative pressurewith height till reaching the maximum negative pressure at the location of thepyramidon at the top of each obelisk. However, the difference in taper anddifference in slenderness between the two obelisks caused the negative pressure onthe Vatican obelisk to be higher than that acting on the Lateran obelisk as theLateran obelisk has a higher taper and is more slender as shown inFigs. 12a and 13a.

On the other hand, the variation in the pressure coefficients in thevertical plane parallel to the wind are very similar for both obelisks showing anincrease in the positive pressure with height in the windward direction. However,the difference in the taper of the pyramidon in each obelisk caused the pressurecoefficients to significantly differ between the two obelisks at the location of thepyramidon as shown in Figs. 12b and13b.

Meanwhile, the variations of the total wind forces with the angles ofprojection for different velocities are plotted for the Vatican obelisk and theLateran obelisk and shown in Figs. 14 and 15 respectively.As expected for both obelisks, the magnitudes of the forces increase with theincrease in wind velocity and the major component of the force for all of the casesis in the x-direction as shown in Figs. 14a and 15a which is thedirection parallel to the wind load itself. When examining the total forces in thez-direction (perpendicular to the wind load) shown in Figs. 14b and 15b, the forces were clearly of lower order than the forces parallelto the wind direction and these transverse forces were 0 for the two symmetric casesof 0o and 45o angles which isexpected to happen. Meanwhile, the general trend for the resultant forces shown inFigs. 14c and 15c show a general trend of the forces increasingwith the angle of projection angle reaching a maximum value at an angle of45o causing this angle to be the most critical loadcase for both obelisks. Also as expected when comparing the resultant forces actingon the Vatican obelisk shown in Fig. 14c to the resultant forces acting on the Lateran obelisk shownin Fig. 15c, the forces acting onthe Lateran obelisk were significantly higher than the forces acting on the Vaticanobelisk which could be attributed to the fact that the Lateran obelisk was larger insize and hence larger in surface area subject to wind load.

Fig. 14
figure 14

Variation of the total wind load acting on the Vaticanobelisk with the angle of projection for different windvelocities

Full size image
Fig. 15
figure 15

Variation of the total wind load acting on the Lateranobelisk with the angle of projection for different windvelocities

Full size image

Furthermore, the variation of the total resultant force with the windvelocity was directly proportional to the square of the velocity as shown inFig. 16. This is considered asproof that the results produced by the CFD are valid and sufficiently accurate interms of obeying the principal relationships between the wind speed and the windforce.

Fig. 16
figure 16

The variation of the total resultant force with differentvelocities

Full size image

The forces produced by the CFD were used in the structural analysis ofeach of the two obelisks when subject to a 30 m/s wind event at the fourdifferent projection angles. The maximum force reached at an angle of45o caused the highest normal stresses within each ofthe two obelisks as shown in Fig. 17where the normal stresses increased with the increase in projection angle.Meanwhile, and despite the fact that the Lateran obelisk had a larger cross-sectionand was expected to have lower stresses, the maximum normal stresses within theLateran obelisk was higher than the maximum normal stresses in the Vatican obelisk.This is due to the higher wind forces acting on the Lateran obelisk compared to theVatican obelisk in addition to the fact that the Lateran obelisk was taller hencethe bending moments experienced by this obelisk were larger than its Vaticancounterpart. However, all of these maximum stresses were even less than0.35 MPa which is significantly less than the strengths reported by[6]. That could explain how thesestructures existed for millennia and survived natural disasters along that longperiod of time as the material used to carve these structures had a strength thatwas higher than any stresses these structures could encounter while loadedlaterally.

Fig. 17
figure 17

The maximum normal stresses within the two studied obelisksdue to different wind load projection angles

Full size image

4 Conclusions and recommendations

In lieu of the results found, the following could be concluded:

  • The CFD models have proven to be valid in terms ofproducing symmetric pressure coefficient contours and symmetricresultant forces for the symmetric load cases of angles of attack of0o and45o.

  • The CFD models have proven to be valid in terms ofproducing resultant forces that are directly proportional to the squareof the reference wind speed.

  • The slenderness and the taper are the main factors thatgovern the variation in pressure coefficient whether along the verticalor horizontal planes.

  • The variation of the pressure distribution at the top ofthe obelisks significantly varied with the dimensions and angles of thepyramidon existing at the top of each obelisk.

  • The largest resultant forces were achieved at an angle of45o for both of the studiedobelisks.

  • The largest normal stresses were achieved at an angle of45o for both of the studiedobelisks.

  • For all of the four different angles of attack included inthe study, the total resultant forces acting on the Lateran obelisk werelarger than the total resultant forces acting on the Vatican obeliskmainly due to the larger surface area of the Lateran obelisk.

  • The structural analysis revealed that for all of the fourdifferent angles of attack included in the study, the normal stresseswere significantly lower than the strength of the Red Aswan granite fromwhich the obelisks were carved.

  • For all of the four different angles of attack included inthe study, the structural analysis revealed that both obelisks couldsafely withstand boundary layer winds with reference speeds as high as30 m/s.

The following recommendations could be drawn from the performedstudy:

  • The responses of the obelisks to wind load need to beexperimentally studied in wind tunnel tests in order to study theeffects of turbulence, presence of surrounding structures and terrainvariations.

  • The responses of the obelisks need to be studied under highintensity winds such as downbursts and tornadoes.

  • Thorough studies need to be performed by engineers incooperation with Egyptologists to further understand how the ancientEgyptians designed such obelisks to withstand such winds further thanhow they designed them.