Ying Chen MSC Report
Ying Chen MSC Report
Ying Chen MSC Report
Ying Chen
Ying Chen
The undersigned hereby certify that they have read and recommend to the European
Wind Energy Master - EWEM for acceptance a thesis entitled Bucket Foundation in
Clay for OWT Subjected to Combined Cyclic Loads by Ying Chen in partial
fulfillment of the requirements for the degree of Master of Science.
Supervisor:
Prof. A. Metrikine of TU Delft
Supervisor:
Prof. H.P. Jostad of NTNU
Daily supervisor:
K.S. Skau of NGI
Reader:
Dr.Ir. R.B.J. Brinkgreve of TU Delft
Reader:
Dr.Ir. W. Broere of TU Delft
Summary
The need for cost reduction of OWT has been a great priority for the wind industry since
the early beginning of the development of OWT. The foundation, account for around one
third of total cost of offshore wind turbine system, is an important part to be optimized
to bring down the cost of OWT systems. Bucket foundations have been given much
attention recently due to its low cost and easy installation method. However, the cyclic
loading conditions experienced by OWT should be investigated more before applying it
to commercial wind farms.
This thesis project investigates the bucket foundation in clay subjected to cyclic loads.
The cyclic irregular loading history in this thesis was represented by the equivalent number
of cycles Neq with a certain level of cyclic to average stress ratio cy /a which can be
derived using cyclic shear strain accumulation procedure developed in NGI. The first
part is the bearing capacity analysis of bucket foundation in clay. The second part is to
produce the displacement contour diagrams.
The bearing capacity analysis used in-house software BIFURC from NGI. The normalized
failure envelopes were investigated for different degree of cyclic degradation expressed by
Neq , soil profiles and bucket foundation geometry h/D. The two soil profiles used in this
thesis project were constant cyclic shear strength with depth for OCR=40 and linearly
increase cyclic shear strength with depth for OCR=1 & 40. Drammen clay was used as
basis for the evaluation.
The mobilized ultimate capacity were plotted for pure horizontal load, pure vertical load
and pure overturning moment. The loads were applied in the decoupling point, meaning
that the pure horizontal load fail under pure sliding. It was found that the three load
component were equally influenced by cyclic degradation and that the increase and re-
duction in the total failure as function of Neq and Fcy /Fa could be defined by one curve
for a given OCR.
The cyclic effect during combined loading was studied through normalized failure en-
velopes. The normalized failure envelopes do not change for changing cyclic loading
history for the same soil profile. The normalized failure envelopes changed in HM load
plane due to different h/D and different soil profiles. The failure envelopes describes a
v
vi Summary
full surface in the 3D load space HVM. Besides, the failure envelopes could be described
by formulas suggested by Gourvenec and Barnett with some modifications for the failure
envelopes in HM plane.
Displacement contour diagrams were produced for a wide range of load combination in
HVM load space. MATLAB was used to interpolate the results from bearing capacity
analysis to produce the displacement contour diagrams. The displacement contours gives
the stiffness for bucket foundation under cyclic loading.
In the end of this thesis, a procedure was made to account cyclic loading to determine the
bearing capacity and stiffness for bucket foundation for OWT. The cyclic loads as input
parameters are represented by Neq and Fcy /Fa . The other input parameters are bucket
foundation geometry h/D and soil profile.
Suggestions for further work can be found in Section 6.2.
Acknowledgements
I would like to express my deep gratitude to Kristoffer Skjolden Skau and Professor Hans
Petter Jostad from NGI, for proposing this thesis topic, for their daily supervision and
for reviewing on my report.
I also would like to appreciate the valuable comments from Professor Andrei Metrikine,
Dr.Ir.Ronald Brinkgreve and Dr.Ir.Wout Broere from TU Delft. I want to thank to all
the EWEM coordination members and 2012-2014 EWEM cohorts for their help during
last two years.
I would like to thank NGI for offering the opportunity for me to join the professional
geotechnical working environment, the chances to attend many lectures offered from
geotechnical experts around the world. I also want to thank all the friendly employ-
ees and their weekly cake from CGM in NGI.
Last but not least, I would like to thank to my parents for their encouragement and
support for my study.
vii
viii Acknowledgements
Nomenclature
Latin Symbols
Fa Average load on the bucket foundation
Fcy Cyclic load on the bucket foundation
Hult Ultimate horizontal capacity
0
K0 Coefficient of earth pressure at rest
Mult Ultimate overturnning moment
Neq Equivalent number of cycles
NF Number of cycles
su Undrained static shear strength
sC
u Undrained static shear strength under triaxial compression
sDSS
u Undrained static shear strength under direct simple shear
sE
u Undrained static shear strength under triaxial extension
ucy Cyclic pore pressure component
up Permanent pore pressure
v0 Poissons ratio
Vult Ultimate vertical capacity
d Wall thickness of bucket foundation
u Horizontal displacement
w Vertical displacement
Greek Symbols
ix
x Nomenclature
Abbreviations
1P First excitation frequency of wind turbines,Rotor rotational frequency
2D 2 Dimensional
3D 3 Dimensional
3P Second excitation frequency of wind turbines,Rotor blade passing frequency
ALS Accidental Limit State
CGM Computational Geomechanics
C Triaxial Compression
DNV Det Norske Veritas
DSS Direct Simple Shear
D Diameter
EU European Union
EWEM European Wind Energy Master
E Triaxial Extension
FEA Finite Element Analysis
Nomenclature xi
Summary v
Acknowledgements vii
Nomenclature ix
xiii
xiv Contents
3 FEA in BIFURC 17
3.1 BIFURC introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
3.1.1 2D simplification with side shear roughness factors . . . . . . . . . 18
3.2 Comparison of BIFURC with PLAXIS . . . . . . . . . . . . . . . . . . . . 19
3.3 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
References 46
C MATLAB code 57
C.1 MATLAB code for batch processing . . . . . . . . . . . . . . . . . . . . . 57
C.2 MATLAB code for plotting 3D failure surface . . . . . . . . . . . . . . . . 60
C.3 MATLAB code for plotting 3D displacement contour diagrams . . . . . . 63
4.1 Three ultimate loading conditions for linear increase cyclic shear strength 21
4.2 Three ultimate loading conditions for constant cyclic shear strength . . . 21
xvii
xviii List of Figures
xxi
xxii List of Tables
Chapter 1
Cyclic loading of soils and foundations has been investigated by many researchers around
the world. This research include different approaches to describe cyclic soil behaviour,
such as cyclic constitutive models, foundation models, pore pressure accumulation di-
agrams etc. The accumulation procedure developed at NGI is one such method. The
method has been widely used and validated through several model tests and field tests.
However, the method is developed on soil element level and there is still need to describe
effects of cyclic loading on foundation level. This thesis investigates the effect of cyclic
loading on the foundation capacity and the foundation stiffness for bucket foundations in
clay. Combined loading of vertical, horizontal and moment loads were considered.
A common way to express foundation capacity for combined loading is to establish failure
envelopes. There is however to the authors knowledge no existing failure envelopes were
adjusted for cyclic loading.
Foundation stiffness can be described using displacement contour diagrams which is sim-
ilar to strain contour diagrams used for soil element.
1.1.2 Limitations
The limitations for bucket foundation in clay subjected to combined cyclic loads investi-
gated in this thesis are as follows:
1
2 Introduction and literature study
- The combined cyclic loads are assumed in-plane loading condition, namely the
three loads are applied to the foundation in the same plane.
- Two idealized soil profiles are used in this thesis. One soil profile has linear
increase cyclic shear strength with depth for both normal consolidated
clay(OCR=1) and heavily over consolidated clay(OCR=40) and another soil
profile has constant cyclic shear strength with depth for heavily over consolidated
clay(OCR=40).
1.1.3 Objectives
The first object is to investigate the cyclic loading effects on pure loading situations
which are pure horizontal load, pure vertical load and pure overturning moment load.
The ultimate cyclic capacity for horizontal, vertical and moment will be compared to the
corresponding ultimate static capacity.
The second objective of this thesis is to produce failure envelopes for bucket foundation
subjected to combined cyclic load for different soil profiles and different bucket foundation
geometries. This thesis aims to investigate whether the failure envelopes will be changed
due to the cyclic loads, the bucket foundation geometry and the soil profile.
The third objective of this thesis was to produce the displacement contour diagrams
for varying cyclic load history and bucket foundation geometry h/D. The displacement
contour diagram gives the soil-structure stiffness which can be used in bucket foundation
design phase.
If possible, it was also a final objective to systemize and organize the results so it could
be used efficiently for approximations of bucket capacity and foundation stiffness of a
foundation subjected to cyclic loading.
The rest of this chapter is the literature study on the bucket foundation for OWT sub-
jected to combined cyclic loads.
Offshore wind power has developed fast the last decade in Europe. Figure 1.1 shows that
the annual offshore wind power installations (red columns) has been increased from 1%
of the annual EU wind energy installation in 2001 up to 14% in 2013. Because of the
total offshore wind power installation targets of 40 GW produced within 2020 for EU, it
is expected to drive the offshore wind energy production to increase in the coming years.
Figure 1.2 shows OWT are installed in continually deeper water. Deeper water will
impose higher loads on the foundation which may change the foundation methods for
OWT. Figure 1.3 shows typical foundation concepts for OWT suitable for different water
depth. The wind industry is forced to reduce the cost of offshore wind to make it a more
competitive renewable energy source. Figure 1.4 gives the breakdown cost for typical
wind turbines and the cost of foundation is around 30%. Reducing foundation cost will
therefore significantly reduce the overall cost for OWT.
1.3 Foundation types for OWT 3
Figure 1.3: Different foundation types for OWT from Byrne [6]
Gravity based, mono-pile and multi-pile foundations have been applied in offshore wind
farms (OWF) while caisson type foundations have not been used in commercial wind
farms until now. However, one prototype of wind turbine with bucket foundation has
been installed in the test field in Frederikshavn in December 2002[1]. Figure 1.5 shows
4 Introduction and literature study
that the dominating foundation type up to now has been mono-piles. Different foundation
concepts are described briefly in the following sections.
80%
70%
60%
50%
40%
30%
20%
10%
0%
Foundation types for offshore wind turbines
Most of the gravity based foundations (Figure 1.3) are made of concrete with or without
skirts which can be ballasted with sand to obtain sufficient dead weight to resist the
horizontal load and overturning moment. It was first applied for OWT in 1990s at the
world first offshore wind farm Vindeby in Denmark.
The gravity based foundation has lots of advantages and is still developing to be suitable
in deeper water and to be able to easy install. The suitable water depth of gravity based
1.3 Foundation types for OWT 5
foundation can vary from several meters to hundred meters. The installation generates
low noise which has small environment impact. It is durable and has low maintenance
requirement. It can be designed to support large wind turbines like 8MW and it has low
fatigue sensitivity.
However, as OWT go to deeper water and harsh environment, the gravity foundation
also faced lots of challenges, such as heavy installation vessels which can transport and
install the giant gravity foundation. The decommissioning also caused troubles for gravity
foundations as it is not easy to take away.
Mono-pile foundation shown in Figure 1.3 is normally made of cylindrical steel where the
tip is connected to the transition piece or to the tower directly and the end is penetrated
into the seabed to provide sufficient capacity. It is a relative simple structure and is the
most common foundation type in OWF. The typical water depth for mono-pile foundation
is from 0 to around 35m. The installation is typically carried out by hammer driving. The
noise generated by the driving activities has raised environment concerns. The fabrication
and installation challenges of large diameter piles required for large water depth may limit
the use in the future.
The mono-pile foundations are flexible and the natural frequency of the OWT should be
designed to avoid the resonance with the rotor rotational frequency 1P and blades passing
frequency 3P as well as wave frequency. This may represent a challenge in deep water
and for large OWT with higher loads.
Multi-pile foundation can be used with jacket and tri-pile support structures. The multi-
pile support structures are generally stiffer and the resonance is less problematic for
multi-pile foundation compared to the mono-pile foundations.
In addition, the loading regime, which is different when compared to mono-pile founda-
tions. The global overturning moment at mudline is taken up by vertical compression and
extension force. The tri-pile foundation has been used in BARD Offshore 1 wind farm in
Germany in 2013.
Bucket foundation is also called suction caisson or skirted foundation. Bucket foundations
are usually a cylindrical segment of steel or concrete. The cylinders are open at the bottom
and closed at the top. It has been used wildly in offshore oil&gas industry in almost 30
years. The applicability of bucket foundations for OWT is still on research level and
not taken into commercial use. Relative low costs and easy installation make bucket
foundation interesting for OWT. Its main characteristic are as follows:
However, it can be very sensitive to the soil conditions compared to other foundation
types. The detailed discussion of bucket foundation for OWT will be given in the following
sections.
According to DNV[9], the geotechnical investigation should include both laboratory test
and in-situ testing. The main geotechnical data should be included as follows:
- Shear strength and deformation properties, as required for the type of analysis to
be carried out;
In terms of soil shear strength, DNV mentions that the shear strength should be deter-
mined properly due to cyclic loading. The cyclic shear strength might be smaller than
static shear strength due to the large cyclic horizontal force and overturning moment[10].
OWT are normally subjected to a combination of wave, current and wind loads. The
loads can be divided into permanent loads and environment loads.
The permanent loads are as follows[9]:
- wind loads;
- hydrodynamic loads induced by waves and current, including drag forces and
inertia forces;
- earthquake loads;
- current-induced loads;
- tidal effects;
1.5 Bucket foundation for OWT 7
- marine growth;
The cyclic environmental loads dominant the loads applied to OWT. The cyclic loads
consist of both wave and wind loads, making the frequency content complex. The cyclic
nature of the loads significantly influence the soil and foundation behaviour.
Bucket foundation shown in Figure 1.6 has been discussed as a foundation method for
OWT for many years[1]. It has been used widely for offshore structures like suction piles
and suction anchors for oil & gas platforms and mooring systems, respectively. Never-
theless, the design of bucket foundation for offshore wind turbine is different compared
to foundation design for conventional offshore structures where the dominate force is
self-weight. The dominant forces for OWT are horizontal force and overturning moment.
Figure 1.6: Bucket foundations for a wind turbine: left: monopod, right: tripod[1]
The critical load applied to a mono bucket foundation and multi-leg support structure
are different. The critical loads are horizontal load and overturning moment for the mono
bucket foundation. While the tension and compression loads may occur at some legs
for multi-leg support structures in addition to horizontal loads and overturning moment.
Figure 1.6 shows load regimes for bucket foundations with monopod support structure
and multi-leg support structure, respectively.
8 Introduction and literature study
The design of bucket foundation can be divided into three phases which are design basis,
conceptual design and detailed design. This project cannot include all of them due to the
time limitation.
This thesis focus on the ultimate limit state(ULS) design for undrained condition which
is part of the conceptual design to determine the dimensions of the bucket foundation.
The other aspects in the conceptual design phase are loading combinations, other limite
states(SLS,ALS,FLS) and the installation study[1].
The bearing capacity in ULS analysis of bucket foundation is usually determined either
by limiting equilibrium methods or FEM[11].
For the undrained condition, the capacity is governed by an undrained cyclic shear
strength, the load combinations and the bucket foundation geometry[11]. The deter-
mination of undrained cyclic shear strength is the first part of the undrained bearing
capacity analysis. Then the undrained bearing capacity can be analysed for a specific
bucket foundation geometry under certain loading combinations.
1.6 Summary
This chapter gives information about this thesis topics and objectives. Besides, the liter-
ature study also presented in this chapter.
The literature study includes the offshore wind energy development and different foun-
dations for fixed OWT. The features for different foundations were included to compare
with the bucket foundation which was investigated in this thesis.
For the bucket foundation, the development of bucket foundation and advantages of using
bucket foundation were enclosed in this chapter. The loading regimes were introduced
here for using bucket foundation for mono-pod and jacket support structure for OWT
respectively.
The next chapter presents the soil behaviour under cyclic loading. A NGI accumulation
procedure method for predicting soil behaviour under cyclic loading was enclosed in the
next chapter.
Chapter 2
Cyclic loading is important for offshore foundation design as it may reduce the soil shear
strength as well as stiffness[12].
The reason for reduced shear strength and stiffness is that the cyclic loading tends to
breakdown the soil structure and increase the pore pressure during undrained conditions.
The left picture in Figure 2.1 shows the stress history with constant average shear stress,a
and constant cyclic shear stress,cy . The right picture in Figure 2.1 gives the pore pressure
change due to the stress history.
According to Figure 2.1, the permanent pore pressure,up increases with time in addition
to a cyclic pore pressure component, ucy fluctuates on top of the permanent pore pressure.
The increase of pore pressure leads to increasing cyclic and average strain.
Figure 2.1: Pore pressure as function of time under cyclic loading. up , permanent pore
pressure,0 , initial consolidation shear stress,ucy , cyclic pore pressure[12]
Soil response can be assumed to be perfectly undrained for clay during a storm period.
Figure 2.2 presents the development of cyclic shear strain and average shear strain for two
different stress ratio. The shear failures are large cyclic shear strain in the left and large
average strain in the right. The shear failure defined by either large average shear strain or
large cyclic shear strain or combination of two, depending on the cyclic to average stress
ratio, cy /a applied to the soil element[13]. NGI uses shear strain level of 15% for either
the average shear strain or the cyclic shear strain as the failure criteria. Figure 2.3 shows
9
10 Soil behaviour under cyclic loads
Figure 2.2: Shear failure, (left) Large cyclic strain, (right) Large average strain
the failure contours which defined by number of cycles, NF , and various combination
of a and cy . The evolution of cyclic strain cy and average strain a are indicated on
Figure 2.3 with the dashed lines.
Figure 2.3: Number of cycles to failure, NF for various combinations of a and cy in triaxial
tests on normally consolidated Drammen Clay[13]
It is easy to use Figure 2.3 to determine cyclic shear strength based on a regular loading
history with constant cyclic stress and constant average stress. However, real load histories
2.1 Cyclic shear strain accumulation procedure for predicting soil behaviour
11
for OWT are irregular and difficult to directly apply to Figure 2.3 to calculate the cyclic
shear strength. Accumulation procedures[2] is therefore developed to allow the cyclic
strength to be calculated for cyclic stress histories with variable amplitude and average
stress.
Before applying the accumulation procedures, the irregular loading histories are organized
into many parcels with constant cyclic shear stress. Then, those parcels are applied to the
strain contour diagrams or pore pressure contour diagrams by assuming the soil element
remember the effects from previous loading parcels. Between different loading parcels, the
shear strain or the pore pressure are assumed constant for shear strain contour diagram
and pore pressure diagrams, respectively. After applying all the load parcels, the final
number of cycles is the equivalent cycles Neq , the final cyclic stress cy and average
shear stress a , which then is assumed to represent the full load history. Applying Neq
together with cyclic to average ratio cy /a in Figure 2.3, the cyclic shear strength can be
determined. The same procedure can be used to calculate the cyclic shear strength for
both triaxial and DSS stress path.
The shear strain accumulation and pore pressure accumulation has been used both for
clay and sand. The shear strain accumulation procedure uses permanent or cyclic shear
strain as the memory parameter to do the accumulation while pore pressure accumulation
procedure uses pore pressure as the memory parameter.The cyclic shear strain accumula-
tion procedure is suitable for clay because shear strain can be measured accurately while
pore pressure accumulation might difficult to precisely measure. The strain contour or
the pore pressure contour diagrams are determined based on cyclic triaxial and DSS lab-
oratory test. The following part illustrates how to use the cyclic shear strain contour
diagram to determine the equivalent number of cycles Neq .
The equivalent number of cycles Neq with a corresponding cyclic and average shear stress
is defined to have the same effect as the actual stress history. Neq is determined by
the cyclic shear strain accumulation procedure developed by Andersen [15]The example
below describes the process of predicting the soil behaviour using cyclic shear strain
accumulation with irregular cyclic loads [2].
1) The irregular load history which is shown in the left of Figure 2.4 can be
transformed into the table in right of Figure 2.4 using rain flow counting method;
2) From Figure 2.5, the cyclic shear stress 36kPa with 100 cycles can be applied in
the cyclic shear strain contour diagram which is shown in path A to B. Then the
12 Soil behaviour under cyclic loads
shear stress increase to 41kPa with constant cyclic shear strain at point C. Point
D is the start point of next load parcel which is 41kPa with 35 cycles and this load
lead to point E which is 70 cycles with cyclic shear strain of 1.5%. This 70 cycles
is the equivalent cyclic numbers Neq . The equivalent cycles Neq = 70 with the
cyclic shear stress in the last load parcel of cy = 41kP a has the same effect as the
irregular load histories shown in Figure 2.4 and can be used later for determining
cyclic shear strength.
The second load parcel start from point D rather than point C was because the
cyclic shear strain increased due to change to the higher cyclic shear stress which
can be determined from the left picture of Figure 2.5. Then the cyclic shear strain
changed from pont C point D without shear strength degradation.
The procedure can be used with both permanent strain and pore pressure as
memory or accumulation parameter. However, a small correction has to be
included when cyclic shear strain as parameter. the correction account for
increasing cyclic shear strain from increasing shear stress without degradation.
The cyclic shear strain changes as cyclic shear stress changes even without cyclic
degradation as shown in Figure 2.6[2].There is an immediate increase in cyclic
A to B .
shear strain cy,i when the cyclic shear stress increases from cy cy
To find Neq at failure, the strain contours have to be scaled down so the accumulation
ends up at failure. This means the material factor account for uncertainties in the full
cyclic behaviour and not only the final maximum load.
Some soil elements along the potential failure surface below bucket foundation experience
different stress conditions which are shown in Figure 2.7. Therefore, the anisotropic cyclic
shear strength determine the failure of bucket foundation. The strain compatibility needs
to be considered when using anisotropic stress dependent undrained shear strength[13].
It means that the failure cyclic shear strain and average shear strain cy + a on the
2.2 Implimentation of soil cyclic behaviour in foundation design 13
A B
Figure 2.6: cyclic shear stress-strain curves when cyclic shear stress is increased from cy to cy
Figure 2.7: simplified stress conditions along potential failure surface in the soil for the
foundation subjected to combined cyclic loading
failure surface should be the same in a different stress element. This will lead to stress
redistribution, and may change the cy /a , along the potential failure surface[13].
The FEA program BIFURC will be used for analysis in this thesis project. BIFURC is
developed in NGI. It is a computer program using FEM for determining of undrained
bearing capacity of embedded structures. BIFURC is a 2D program that takes into
account 3D effects using side shear roughness factors[11].
Drammen Clay contour diagrams has been used in study because of its extensive database
of cyclic tests[12, 14]. Cyclic soil shear strength might be lower or higher than static
shear strength due to cyclic loading. The cyclic shear strength was determined prior to
the bearing capacity analysis. The following section describes the determination of cyclic
shear strength.
Figure 2.8 shows the compression, extension and DSS static shear strength for Drammen
Clay with different OCR. The diagram can be approximated by the SHANSEP formula[16]
as follows:
0
su /vc = A OCRB (2.1)
0
Where, vc is the vertical effective consolidation stress; A, B are the coefficients need to
be decided; OCR is over-consolidation ratio.
Coefficient A is obtained from Figure 2.8 by setting OCR=1. The B coefficient was
determined from the average results giving different OCR. The normalized anisotropic
strength is shown in Table 2.1 as function of OCR.
The strength and stress-strain curves were derived using the cyclic test data on Drammen
Clay which has been well documented from K.H. Andersen[14]. Simplified procedures were
used in this thesis project to derive the cyclic shear strength and stress-strain curves.
The cyclic shear strength was defined as, f,cy = a,f + cy,f , where a,f is average shear
strength at failure, cy,f is the cyclic shear strength at failure. The combination of average
shear strength and cyclic shear strength is given in Figure 2.9 for different number of cycles
to failure NF and different cyclic to average stress ratios cy /a . Cyclic shear strength
and average shear strength at failure can be obtained for specific Neq and cy /a .
2.3 Soil cyclic shear strength and stress-strain curve 15
Figure 2.8: Undrained DSS and triaxial monotonic shear strength of Drammen Clay as
function of OCR[10]
The cyclic shear strength was derived assuming that the strain from the given cy /a ratios
were governing for the problem. The average shear strain at failure taua,cy and cyclic shear
strain at failure cy,f which is obtained from DSS contour diagram (Figure 2.9) based on
the cy /a ratio under consideration were used in the compression and extension contour
C = C + C
diagrams(Figure 2.9) to acquire compression cyclic shear strength f,cy a,f cy,f and
E E E
extension cyclic shear strength f,cy = a,f + cy,f . The normalized cyclic shear strength
Figure 2.9: Results from cyclic DSS test and cyclic triaxial test on Drammen Clay with
OCR=4[10]
for all the cases considered in this thesis project are given in Appendix A.
The stress strain curves were generated to represent the characteristic stress paths in the
contour diagrams. The path should at the 0 relevant.
Figure 2.10 shows one of the stress-strain curves used in this thesis project. All the
stress-strain curves used are shown in Appendix B.
16 Soil behaviour under cyclic loads
D S S
C o m p r e s s io n 1 .5
E x te n s io n
1 .0
D S S
u
o r c y / s
0 .5
c
u
c y / s
0 .0
-2 0 -1 5 -1 0 -5 0 5 1 0 1 5
-0 .5
-1 .0
[% ]
2.4 Summary
This chapter shows the theory to predict soil behavior under cyclic loading developed in
NGI which is used in this thesis. This cyclic shear strain accumulation procedure changes
the irregular loading history into an equivalent number of cycles Neq with a certain level
of cyclic shear stress cy and average shear stress a . It means that an irregular loading
history with many different loading amplitude and cycles can be represented by one Neq
with a specific cy /a using the cyclic shear strain accumulation procedure.
The next chapter introduces the in-house software BIFURC developed in NGI to do the
bearing capacity analysis for general cyclic loading conditions which was represented by
varying Neq and cy /a .
Chapter 3
FEA in BIFURC
The bearing capacity of the buckets were computed with BIFURC. BIFURC is a 2D
program using FEM by taking 3D effects into consideration using side shear roughness
factors. BIFURC can produce the 3D bearing capacity envelop efficiently relative to other
full 3D programs.
BIFURC is a finite element program for failure analysis of materials with anisotropic
shear strength as well as consolidation analysis[3].
The basis of the material model used are as follows:
- An unloading modulus equal to the initial shear modulus Plane strain condition[3].
The soil elements used in BIFURC are 8-node isoparametric Serendipity element with
numerical integration of 2 2 Gauss points[3]. This element is different than the normal
plane element since it has thickness in the y-direction and shear stresses on the two
surfaces perpendicular to the y-direction (Figure 3.1). The interface element is 6-node
isoperimetric element (Figure 3.2) and it is used to model sliding between two 8-node
elements.
17
18 FEA in BIFURC
PLAXIS is developed from PLAXIS BV in The Netherlands. And all the relevant infor-
mation are well documented in their manuals which can be download from their home
website which is http://www.plaxis.nl/.
BIFURC bring 3D effects into calculation using side shear roughness factors which are
shown in Figure 3.3. There are three side shear roughness factors in BIFURC which are
reduction factor for shear on structure outside area bs , reduction factor for shear on
structure inside side area, reduction factor for shear on soil side areas ss . bs is at the
two plane vertical sides of the bucket foundation and ss is at sides of the active and
passive zones. Figure 3.3 shows the side shear roughness factors at the two plane sides of
the failure body for shallow foundation.
3.2 Comparison of BIFURC with PLAXIS 19
The side shear roughness factors has been extensively be verified from K.H. Andersen and
H.P. Jostad[11]. This section compared the ultimate capacity results using BIFURC with
roughness factors recommended from K.H. Andersen and H.P. Jostad[11] to PLAXIS.
The finite element models were established in PLAXIS 2D, PLAXIS 3D and BIFURC. The
soil profile used in this section was OCR=1. The load condition was a /cy = Fcy /Fa = 1
and equivalent cycles Neq =1. BIFURC calculates the ultimate bearing capacity with side
shear roughness factors shown in Table3.1.
Table 3.1: Side shear roughness factors in BIFURC
The mesh in PLAXIS 2D and PLAXIS 3D are sufficient to produce the ultimate bearing
capacity analysis. The vertical ultimate bearing capacity using BIFURC was compared
to the result from PLAXIS 2D shown in Figure 3.4. The bearing capacity using BIFURC
was 4% higher than the bearing capacity from PLAXIS 2D.
For the moment bearing capacity, the overshoot from PLAXIS 3D compare to BIFURC
was 5% which is shown in Figure 3.5. The bearing capacity using BIFURC with side
shear roughness factors shown in Table 3.1 agreed well with PLAXIS 2D and PLAXIS
3D. The overshoot around 5% can be neglected and may be due to the different element
and interface in BIFURC and PLAXIS.
The agreement was considered to be sufficient for this study and the roughness factors in
Table 3.1 was used in the analysis presented herein.
3.3 Summary
This chapter introduces the material model in BIFURC and elements used in BIFURC.
This chapter also describes how the 2D element can take 3D effects by side shear roughness
20 FEA in BIFURC
3 0
2 5
2 0
(M N )
1 5
u lt
V
1 0
P L A X IS 2 D
B I F U R C w i t h b s = 0 . 5 , s s = 0 . 5
5
0
0 .0 0 .2 0 .4 0 .6 0 .8 1 .0
w (m )
6 0 0 0 0
5 0 0 0 0
4 0 0 0 0
(k N m )
3 0 0 0 0
u lt
2 0 0 0 0
M
B I F U R C w i t h b s = 0 . 5 , s s = 0 . 5
1 0 0 0 0 P L A X IS 3 D
0
0 .0 0 .2 0 .4 0 .6 0 .8 1 .0
u (m )
This chapter presents the normalized bearing capacity envelopes in three load dimen-
sions,VHM. Vult , Hult and Mult , denote pure vertical load, pure horizontal load at ref-
erence point, RP, and pure moment load(Figure 4.1). The RP is the determined as the
point where a horizontal load can be applied and cause a pure horizontal sliding without
rotation. The RP was change from 2/3h for linear increasing f,cy with depth to 1/2h for
constant f,cy with depth shown in Figure 4.2.
Figure 4.1: Three ultimate loading conditions for linear increase cyclic shear strength with
depth
Figure 4.2: Three ultimate loading conditions for constant cyclic shear strength with depth
For combined loading conditions shown in Figure 4.3, horizontal, vertical load and over-
turning moment were applied on the RP.
The analyzed cases in this thesis project are listed in Table 4.1.The soil-skirt roughness
factor is set to 0.65 for all cases.
21
22 Failure envelopes determination
Figure 4.3: General loading conditions for linear increasing f,cy with depth(left) and constant
f,cy with depth(right)
The bucket foundation geometry of embedment depth to diameter ratios h/D are 0.5 and
1. For OWT it is more realistic to use shallow foundations for monopod structure since
the dominant load usually is overturning moment. The model idealization used embedded
solid structure. The failure surface is outside the bucket where the bucket was assumed
to have inside stiffener which prevents the failure mode to go up between the skirts.
The undrained cyclic shear strength are taken from previous section where it considered
loading history and OCR. The undrained shear strength is assumed constant from mudline
down to 1m and linearly increases with depth. Figure 4.4 shows the bucket foundation
geometry and the corresponding assumed soil condition.
In order to plot the full failure envelopes, different load combinations were applied to
calculate the failure load for each load history and OCR. Many load combinations were
considered, and batch processing using MATLAB script were developed to reduce the
calculation time.It writes loads into BIFURC loads file and then activates the executive
file from BIFURC. The user only need to input predefined approximate loads for pure
horizontal load, pure vertical and pure overturning moment. MATLAB writes the ulti-
mate loads into the excel file. The rest general load combinations were based on the three
ultimate loads. MATLAB reads again the load combination excel file and write each load
combination into load file for BIFURC and call BIFURC to calculate.
The results were written into the workspace in MATLAB. Finally, all the load-displacement
results were written into one of the sheets in an excel file after calculation finished for all
the load combinations for one soil profile. The MATLAB code for the batch processing
is shown in Appendix C.1.
The effects of cyclic loads on the ultimate capacity were investigated by comparing results
from different Fcy /Fa and Neq .The ratios between cyclic ultimate capacity and static
ultimate capacity which are Mult /Mstatic , Hult /Hstatic and Vult /Vstatic are almost the
same for horizontal, vertical and overturning moment loading. The effect of cyclic loading
expressed as Neq and Fcy /Fa can therefore be reduced to two figures shown in Figure 4.5.
This figure also represent the effects on Mult /Mstatic and Vult /Vstatic .
1 .6 1 .6
F a
= 0
F /F = 0 .5 1 .4
1 .4 a c y
F /F = 1 1 .2
a c y
1 .2 1 .0
s ta tic
s ta tic
0 .8
1 .0
/H
/H
u lt
u lt
0 .6
H
N = 1
H
0 .8 0 .4 e q
N e q
= 1 0
0 .6 0 .2
N e q
= 1 0 0
0 .0
1 1 0 1 0 0 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0
N e q
F a
/F c y
4
x 10 1
6
0.9
5
0.8
H [KN] or M [KNm]
4 0.7
H/Hult or M/Mult
3 0.6
0.5
2
0.4
1
0.3
VH plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0
0 0.2
0 0.5 1 1.5 2 2.5 3 3.5 VM plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0
V [KN] or H [KN] x 10
4
0.1 HM plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0
VH plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0
VM plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
HM plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0 V/Vult or H/Hult
Figure 4.6: Failure envelopes for different load planes for h/D=1
Failure envelopes were computed for different OCR, Neq , Fcy /Fa and the bucket geometry
h/D. The input cyclic shear strength can be found in Appendix A. The combined bearing
capacity was normalized by the corresponding ultimate bearing capacity. The shape of
the failure envelopes were first investigated in two load planes. Then the combined three
loading conditions was further applied to the bucket foundation.
This section investigate whether cyclic loading expressed as cy /a and Neq , has influence
the failure envelopes.
The failure envelopes are shown in the left of Figure 4.6 in real number. The normalized
failure envelopes for different Neq and cy /a were plotted in the right of Figure 4.6 in
HM, VH and VM load plane for Neq =1, 10,100 and cy /a = 0, 1, 2,a = 0 respectively.
The normalized failure envelopes are almost identical in the same load plane for different
Neq and cy /a . The cyclic loads do not influence too much on the normalized failure
envelopes for the same bucket foundation geometry h/D.
Figure 4.7 shows the failure envelopes in VH,VM and HM plane not change due to higher
OCR which gives a higher cyclic shear strength gradient.
4.3 Failure envelopes for linear increasing cyclic shear strength with depth25
0.9
0.8
0.7
0.5
0.4 VH plane,OCR=1
VM plane, OCR=1
0.3
HM plane, OCR=1
VH plane, OCR=40
0.2
VM plane, OCR=40
0.1 HM plane, OCR=40
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult or H/Hult
Both failure envelopes and the normalized failure envelopes for h/D=0.5 and 1 were
analyzed in this section for OCR=1, cy /a = 2. Figure 4.8 shows the failure envelopes
in real numbers and normalized failure envelopes. The shapes of the failure envelopes in
VH and VM planes in Figure 4.8 almost identical which might be due to the geometry
change is small. The failure mechanisms in VH plane and VM plane stays almost the
same for bucket foundation geometry h/D=0.5 and 1.
The failure shape in HM plane significantly changes between h/D=0.5 and 1 which can
be seen from Figure 4.8. The coupling between horizontal load and moment load becomes
less when h/D goes small namely more shallow foundation.
5
x 10 1
3
0.9
VH plane, OCR=1
2.5 VM plane, OCR=1 0.8
HM plane, OCR=1
VH plane, OCR=40 0.7
2 VM plane, OCR=40
H [KN] or M [KNm]
H/Hult or M/Mult
HM plane, OCR=40 0.6
1.5 0.5
0 0
0 1 2 3 4 5 6 7 8 9 10 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V [KN] or H [KN] 4 V/Vult or H/Hult
x 10
For OCR=40, the high gradient of cyclic shear strength used in previous section was
somewhat unrealistic and a constant shear strength profile with compression cyclic shear
c
strength f,cy = 100kP a was used herein. As the failure envelopes only change due to h/D
from analysis in previous sections, this section compared the failure envelopes between
constant f,cy with depth and linear increase f,cy with depth and failure envelopes for
different h/D for constant f,cy with depth.
4.4.1 Failure envelopes between linear increase f,cy and constant f,cy
with depth
The failure envelopes for constant cyclic shear strength was compared to linear increase
cyclic shear strength are shown in Figure 4.9. The change between VH and VM plane are
relatively small while the change in HM plane are very large because the scaling factors
on H axes and M axes are very different for linear increase f,cy with depth and constant
f,cy with depth. The trend of failure envelopes seem very similar when looking at the
real number in left of Figure 4.9.
The failure envelopes for h/D=0.5 and 1 were presented here to confirm the same effect
from geometry change of bucket foundation. Figure 4.10 shows the failure envelopes in 3
different load plane for h/D=0.5 and 1. The results agree well with failure envelopes for
linear increasing cyclic shear strength in section 4.3.3.
4.5 Failure envelopes comparison 27
5
x 10 1
3
0.9
VH plane, h/D=1
2.5 VM plane, h/D=1 0.8
HM plane, h/D=1
VH plane, h/D=0.5 0.7
2 VM plane, h/D=0.5
H [KN] or M [KNm]
H/Hult or M/Mult
HM plane, h/D=0.5 0.6
1.5 0.5
0 0
0 1 2 3 4 5 6 7 8 9 10 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V [KN] or H [KN] 4 V/Vult or H/Hult
x 10
The failure envelopes are also plotted using the formula from Gourvenec and Barnett[17].
In their paper, the proposed envelope formula are shown below:
For HVM plane:
h m hm
(
) +( ) + 2( ) = 1 (4.1)
h m h m
For only HV plane:
D D
= 1.30 + 1.05( ) 0.55( )2 (4.4)
B B
D D
) + 0.67( )2 = 0.15 1.45(
(4.5)
B B
Where, D is the embedment depth of the bucket foundation; B is the width of the bucket
foundation The failure envelopes were compared for Neq = 0, 10, 100, cy /a = 0, 1, 2 &
a = 0 for OCR=1 and h/D=1. The results are ouside the failure envelopes in VH and
VM plane from Gourvenec and Barnett which are shown in Figure 4.11. This is because
the failure envelopes from Gourvenec and Barnett was a conservative solution.
The comparison in HM load plane cannot do directly because the reference points are
different from Gourvenec and Barnett which they chose the tip of the bucket foundation
28 Failure envelopes determination
1 1
0.9 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
Gourvenec,S. & Barnett,S.(2011) Gourvenec,S. & Barnett,S.(2011)
0.1 0.1
VH plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0 VM plane h/D=1,OCR=1,Neq=1,10,100, cy/a=0,1,2 & a=0
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
as RP instead of the approximate decoupling point in this thesis project. This comparison
was carried out by transferring the load from decoupling point to the tip of the bucket
foundation. Figure 4.12 shows that the failure envelope from Gourvenec and Barnett
overestimates the capacity in HM load plane. And the blue line was a proposed solution
to shrink the numerical failure envelope around the failure envelopes from this project.
The improved formulas are in the following:
h m hm
(
) +( ) + 3.2( ) = 1 (4.6)
h m h m
D D
= 0.75 + 1.5( ) 0.3( )2 (4.9)
B B
D D
) + 0.9( )2 = 0.2 1.47((4.10)
B B
The failure envelope from Gourvenec and Barnett did not considered the physical sit-
uation that the failure mechanism is dominate by the horizontal movement failure until
overturning moment is big than a certain level.
The 3D failure surface can be plotted using the bearing capacity analysis results under
general loading conditions. Before plotting the 3D failure surface, it is necessary to check
the failure envelopes symmetrical characteristic in each load plane to be able to properly
produce the 3D failure surface.
4.6 Full failure envelopes and 3D failure surface 29
1.5
M/Mult
0.5
h/D=1,OCR=1,Neq=10,cy/a=1
Gourvenec,S. & Barnett,S.(2011)
Proposed improved solution
0
1 0.5 0 0.5 1 1.5
H/Hult
This section analyzes the full failure envelopes in HM,VH,VM load planes.
VH and VM load plane symmetrically in four quadrants because the failure mechanisms
in VH and VM load plane do not change for different direction of the loads. It means
that the failure envelopes in previous section for VH and VM load plane are enough to
express the full failure envelopes.
For HM load plane, the failure mechanism is different between the first quadrant and the
second quadrant. Figure 4.13 shows that the failure mechanism in the first quadrant is
clockwise rotation with right horizontal displacement while in the second quadrant, the
failure mechanisms is clockwise rotation with left horizontal displacement. Figure 4.14
shows the full HM failure envelope for OCR=1,Neq = 20, cy /a = 1.
Figure 4.15 shows the comparison of the full failure envelope in HM load plane and the
main difference is in the first quadrant.
30 Failure envelopes determination
1 .0
0 .8
0 .6
u lt
0 .4
M /M
0 .2
0 .0
-1 0 1
H /H u lt
1 .0
0 .8
0 .6
u lt
M /M
0 .4
0 .2
h /D = 0 .5
h /D = 1
0 .0
-1 .0 -0 .5 0 .0 0 .5 1 .0
H /H u lt
Normally three combined loading situations are experienced by the OWT. The full 3D
failure surface from combined cyclic loading conditions are presented in this section. The
results were produced by a series of analysis run by batch processing script in MATLAB to
calculate the bearing capacity for bucket foundation. The 3D failure surface and several
cross section are presented in this section. The MATLAB script used linear interpolation
function from the output data to derive cross sections of failure surface.The MATLAB
script is in Appendix C.2
4.6 Full failure envelopes and 3D failure surface 31
Figure 4.16 shows cross section of 3D failure surface for a given normalised moment load
in HV load plane. Blue failure envelopes are the cross sections for M/Mult from 0.7 down
to 0.3 with a interval of 0.1. The red failure envelopes is in HV plane with overturning
moment equals zero. Figure 4.17 shows the 3D failure surface.
1 M/Mult=0
M/Mult=1
0.9
M/Mult=0.3 0.7
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1
0
0
1 0.8 0.6 0.4 0.5
0.2 0
1
V/Vult
H/Hult
Figure 4.16: Overturning moment cross section of failure envelopes for OCR=1,Neq = 100,
cy /a = 1
4.7 Conclusions
The bearing capacity under cyclic loading were compared to static capacity. It shows that
the cyclic bearing capacity can be higher or lower than the static capacity which depends
on the loading types and cyclic to average ratio cy /a .
The normalized failure envelopes by giving different Neq and cy /a , bucket foundation
geometry h/D and soil shear strength profiles. The failure envelopes does not change duo
to cyclic loading history. However, it changes for different h/D and soil shear strength
profiles.
The next chapter uses MATLAB to interpolate the bearing capacity analysis results to
plot the displacement contours. The displacement contour diagrams will be produced for
two load space and HVM three load space.
Chapter 5
This chapter interpolates the load displacement results from bearing capacity analysis
in last chapter to plot the displacement contours. The load displacement results are
written into MATLAB workspace. The displacement contour diagrams are produced
using MATLAB code. This code uses linear interpolation the corresponding loads for a
given displacement level.
The displacement is the total displacement including cyclic and average displacement.
The secant stiffness can be calculated from the displacement contour diagrams for a
given load level.This secant stiffness for example can be used as boundary conditions for
buckets on a jacket structure in analysis performed to calculate the load distribution in
the structure.
The failure displacement are 1m for pure horizontal and pure vertical loading and 0.14
rad for moment loading. The load-displacement curve goes flat at 1m displacement and
0.14 rad.
The displacement contours are plotted in VM, VH and HM load plane respectively. Fig-
ure 5.1 shows the displacement contours in VM plane with horizontal load equals zero.
The dash lines with circles shows the rotational radian, where the solid lines with asterisk
shows the vertical displacement.
Figure 5.2 shows the displacement contours in VH plane with overturning moment equals
zero. The horizontal displacement u and vertical displacement w were chose from 0.01m
33
34 Displacement contour diagrams determination
H=0
1
= 0.0010.14 rad
0.9 w = 0.011 m
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult
to 1m with an interval of 0.01m. The most outside line where two lines overlap each other
is the defined failure line in VH plane.
M=0
1
u = 0.011 m
0.9 w = 0.011 m
0.8
0.7
0.6
H/Hult
0.5
0.4
0.3
0.2
0.1
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult
Figure 5.3 shows the displacement contours in HM plane with vertical load equals zero.
The contour expand when the horizontal displacement and rotation become small. One
of the reasons might because the reference point was changed in the reality while it is
5.1 2D dispalcement contour diagrams 35
V=0
1
u = 0.011 m
0.9 = 0.0010.14 rad
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
The displacement contour diagrams were not produced for OCR=40 with linear increasing
shear strength since it is somewhat unrealistic. This section compared the displacement
contours for OCR=1 with linear increasing cyclic shear strength with depth to OCR=40
with constant cyclic shear strength with depth. The displacement contours changed due
to change of OCR. Figure 5.4 to Figure 5.6 show the displacement contours in VM, VH
and HM load planes for OCR=1 & 40, Neq = 10 and cy /a = 1. The displacement
contours changed due to the change of the stress strain curve for different OCR.
The increase of embedment depth to diameter ratio, h/D causes the failure envelopes
change in HM load plane. According to Figure 5.7 to Figure 5.9, the geometry change
also causes the displacement contours change in horizontal displacement and rotational
radians. The increased embedment depth do not influence the vertical displacement
contours because the vertical displacement contours related to the diameter of the bucket
foundation.
36 Displacement contour diagrams determination
H=0 H=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
M/Mult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
= 0.0010.14 rad
w = 0.011 m
0.1 0.1
Failure envelope
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
M=0 M=0
1 1
u = 0.011 m
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
H/Hult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
u = 0.011 m
w = 0.011 m
0.1 0.1
Failure envelope
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
V=0 V=0
1 1
0.9 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
M/Mult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
u = 0.011 m
0.1 u = 0.011 m 0.1 = 0.0010.14 rad
= 0.0010.14 rad Failure envelope
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult H/Hult
H=0 H=0
1 1
= 0.0010.14 rad = 0.0010.14 rad
0.9 w = 0.011 m 0.9 w = 0.011 m
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
M/Mult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
M=0 M=0
1 1
u = 0.011 m u = 0.011 m
0.9 w = 0.011 m 0.9 w = 0.011 m
0.8 0.8
0.7 0.7
0.6 0.6
H/Hult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
V=0 V=0
1 1
0.9 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
M/Mult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
The 3D displacement contours on HV load plane with different constant moment loads
are presented in this section. Figure 5.10 shows the displacement contours in HV load
plane with a constant moment load M = 0.5Mult . The rotation is in a range from 0.001
rad to 0.14 rad with an interval of 0.001 rad. The interpolation method used herein was
linear interpolation function from MATLAB. The corresponding MATLAB code was in
Appendix C.3. The core interpolations in this section has two steps.
1) Interpolate to acquire horizontal load and vertical load for all the combined
loading conditions. In addition, the capacity should also add to previous
parameters which need to do the moment interpolation in next step. The
parameters needed for each load combinations is (,Mload , Hload , Vload ,Mcapa ,
Hcapa , Vcapa );
2) For the same , the parameter sets were divided into 5 groups with 7 points in
each group. For each group, horizontal load and vertical load were calculated by
interpolate moment capacity Mcapa . The parameter sets changed into (,Mload ,
Hload , Vload ,0.5, Hcapa , Vcapa ). Each group generates one point on the displacement
contours for a specific and M/Mult . 5 groups give 5 points on the contour
diagrams. The points on the axis were derived in two load plan which are HM
load plane and VM load plane used the method in plotting 2D displacement
contours diagram from Section 5.1.
The displacement contours for horizontal displacement and vertical displacement were
used the same method. Figure 5.10 shows the displacement contours for ,u and w in
HVM load plane.
The load-displacement curves for a bucket foundation can be determined from Figure 5.10
for a given combined cyclic loading. The detailed procedure for calculating the bearing
capacity and stiffness of bucket foundation was presented in next section.
From previous analysis, the failure envelopes and corresponding displacement contours
obtained can be used for the preliminary design of bucket foundation for OWT subjected
to cyclic loads. This section gives detailed procedure to use the failure envelopes and
displacement contour diagrams which were produced in this thesis report to calculate the
bearing capacity and soil structure stiffness by giving a specific irregular load history.This
procedure is shown in Figure 5.11.
Irregular load history transform into load parcels where each load parcel has a constant
average and cyclic load. Load parcels uses NGI cyclic shear strain accumulation procedure
to get a cyclic and average stress cy /a and equivalent cycles Neq . Those are the input
parameters for the procedure introduced in this section.
5.4 Summary 39
M/Mult=0.5 M/Mult=0.5
1 1
0.9 0.9
=0.001 0.14 rad u = 0.01 1 m
0.8 0.8
0.7 0.7
0.6 0.6
V/Vult
V/Vult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult H/Hult
M/Mult=0.5
1
0.9
w = 0.01 1 m
0.8
0.7
0.6
V/Vult
0.5
0.4
0.3
0.2
0.1
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
The static capacity analysis for pure horizontal HStatic , pure vertical VStatic and pure
overturning moment MStatic should performe at first. The cyclic ultimate capacity for
pure horizontal Hult , vertical Vult and overturning moment Mult can be determined from
Neq and Fcy /Fa . This calculation can be carried out using Figure 4.5.
Therefore, the failure envelopes were determined using Figure 5.12.
Finally, the soil structure stiffness were determined from Figure 5.10.
5.4 Summary
This chapter presented the displacement contour diagrams for different Neq and cy /a ,
h/D and soil profiles. In the end of this chapter, a procedure was introduced to guide
how to use the data produced and plotted in this thesis report.
In the next chapter, the conclusions and suggestions will be presented.
40 Displacement contour diagrams determination
Input parameters:
Neq and Fcy/Fa, OCR, h/D
Perform analysis for
Mstatic, Hstatic, Vstatic
Stiffness - displacement
contour diagrams
Load distribution
Figure 5.11: Procedure for determining cyclic bearing capacity and stiffness for bucket
foundation
0.9
0.8
0.7
0.6
V/Vult []
0.5
0.4
0.3
0.2 M/Mult=0
M/Mult=1
0.1
M/Mult=0.3 0.7
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult []
6.1 Conclusions
The research presented in this report has studied the response of bucket foundation for
OWT subjected to combined cyclic loads. Normal consolidated clay and heavily consol-
idated clay have been considered for bucket foundation geometry h/D=0.5 and h/D=1.
The research has successfully revealed several characteristics of the behaviour of a cycli-
cally loaded bucket foundation. Based on the results, simplified design diagrams showing
foundation failure and foundation displacements are established.
The ultimate cyclic bearing capacity Mult , Hult and Vult were investigated for different
number of cycles Neq and average to cyclic load ratio Fcy /Fa . The ultimate cyclic bearing
capacity can be higher or lower than the ultimate static bearing capacity depending on
Neq and Fcy /Fa . The ultimate cyclic bearing capacity normalized by the ultimate static
corresponding capacity are almost identical for Mult /Mstatic , Hult /Hstatic and Vult /Vstatic
and therefore only one diagram is needed for determining ultimate cyclic bearing capacity
for the same soil profile.
The effects from different parameters on the normalized failure envelopes in VM,VH and
HM plane are studied.
The equivalent number of cycles Neq and cyclic to average load ratio Fcy /Fa do not influ-
ence the normalized failure envelopes significantly. The soil profile with constant cyclic
shear strength with depth influences the failure envelopes in HM plane when compared
to the soil profile with linear increasing cyclic shear strength with depth.
The bucket foundation geometry also effects the failure envelopes in HM plane both for
normal consolidated clay and over consolidated clay.
Since the failure envelopes only depends on the soil profile and bucket geometry, only
four failure envelopes are sufficient to describe the variation considered in this study.
41
42 Conclusions and suggestions
However, this conclusion is only validate for the the foundation studied in this thesis.
The normalized failure envelopes for deep bucket may be more sensitive to increasing Neq
than the findings in this study.
Based on the conclusions above, it is possible to estimate the capacity for any load com-
binations if the static capacity Mstatic , Hstatic and Vstatic , degradation expressed as Neq
and Fcy /Fa are known. The procedure can be used for preliminary design for the range
of parameters considered in this report. Layered soil and other bucket geometries may
give other relations.
The failure envelopes are symmetry in first quadrant in VM and VH plane and symmetry
around first and second quadrant in HM plane which has been discussed in Section 4.6.1.
The coupling of three loads, which are horizontal, vertical and moment loads, in each
load plane are different which can be seen from the failure envelopes in VM, VH and
HM plane. The coupling of the failure load in VH and VM plane are almost the same
while the coupling of failure loads in HM plane is different and significantly influenced by
bucket foundation geometry h/D.
The failure envelopes in VM and VH plane have a good agreement to the proposed
analytical failure envelopes from Gourvenec and Barnett[17].
Displacement contour diagrams have been developed by the similar procedure and based
on Drammen clay stress train data. The trends of the displacement contour diagrams have
a good agreement to the displacement contour diagrams for a spud can from Jostad[18].
The general displacement contour diagrams also fits well in terms of the strain contour
diagrams in soil element level from Andersen[14].
The displacement data are represented in diagrams as total numbers and not normalized,
but can be used for preliminary design of bucket geometries close to the geometries con-
sidered. The load data in the diagrams are however normalized and based on the study of
the normalized failure envelopes is reasonable to assume that the displacement contours
can be scaled according to the same rules by the ultimate pure capacities.
The displacement contour diagrams are sensitive to the cyclic parameters Neq and Fcy /Fa
as well as bucket geometry h/D and soil profiles.
The decoupling between horizontal displacement and moment rotation for h/D=0.5 is
more obvious than h/D=1. The coupling of the displacement contours in VH and VM
plane are more or less the same.
Due to time limitation, this thesis project was not cover all the aspects related to bucket
foundation design subjected to combined cyclic loads for OWT. Therefore, several sug-
gestions for the further work are listed in this section.
The geometry of bucket foundation of h/D =0.5 and 1 were analysed from previous
chapters. The h/D ratio has a significant influence on the normalized failure envelopes and
corresponding displacement contours which means it effects the failure mode and stiffness
6.2 Suggestions for further work 43
of the soil-structure system. Future work should study different bucket geometries and
produce the capacity envelopes and displacement contours.
The displacements given in contour diagrams are the total displacement and not nor-
malized by any parameters. In future work, the displacement contour diagrams may
use a similar normalized approach as used for establishing the failure envelopes. The
displacement may be normalized by e.g. geometry parameters h or D. This will make
the displacement contour diagrams more general and applicable to more situations. To
incease the possible applicability, it is also necessary to produce displacement contours
for other cyclic to average load ratios Fcy /Fa and OCR values between 1 and 40.
The failure envelopes and displacement contours can be analysed in the future for bucket
foundation in layered soil.
The soil-skirt roughness factors were the same which is =0.65 for all the cases studied
in this thesis. A different parameter is not expected to violate the main conclusions but
the effect should be studied and qualified.
44 Conclusions and suggestions
References
[1] Lars Bo Ibsen, S Liingaard, and Sren A Nielsen. Bucket foundation, a status.
Proceedings of the Copenhagen Offshore Wind, 2005.
[2] Knut H Andersen and R Dyvik. Clay behaviour under irregular cyclic loading. 1992.
[3] H.P. Jostad. Bifurc - version 3 undrained capacity analyses of clay. Technical report,
Norwegian Geotechnical Institute, 1997.
[4] European Wind Energy Association et al. Wind in power, 2012 european statistics,
2013. Available online a t http://www. ewea. org/fileadmin/files/library/publication-
s/statistics/Wind in power annual statistics 2012. pdf, 2012.
[5] Crown Estate. Offshore wind cost reduction: Pathways study. Study in partnership
with the Secretary of State for Energy and Climate Change in the UK, 2012.
[6] Byron Byrne. Foundation design for offshore wind turbines. Technical report, Uni-
versity of Oxford, 2011.
[7] IRENA. Renewable energy technologies: Cost analysis series. Technical report,
International Renewable Energy Agency, 2012.
[8] Lars Bo Ibsen and Rune Brincker. Design of a new foundation for offshore wind tur-
bines. In Proceedings of The 22nd International Modal Analysis Conference (IMAC),
Detroit, Michigan, 2004.
[9] Det Norske Veritas. Design of offshore wind turbine structures. Det Norske Veritas,
2013.
[10] KH Andersen. Cyclic clay data for foundation design of structures subjected to
wave loading. In Proceedings of the International Conference on Cyclic Behaviour
of Soils and Liquefaction Phenomena, CBS04, Bochum, Germany, volume 31, pages
371387, 2004.
45
46 REFERENCES
[11] Knut H Andersen, Hans Petter Jostad, et al. Foundation design of skirted founda-
tions and anchors in clay. In Offshore Technology Conference. Offshore Technology
Conference, 1999.
[12] Knut H Andersen. Bearing capacity under cyclic loading-offshore, along the coast,
and on land. the 21st bjerrum lecture presented in oslo, 23 november 2007 this paper
represents the written version of the 21st bjerrum lecture. while it has been edited for
the present publication, it retains the general structure of the original lecture, which
was intended for a general geotechnical audience. the bjerrum lecture is presented
in oslo in alternate years by the norwegian geotechnical society with the support of
the bjerrum memorial fund (laurits bjerrums minnefond). Canadian Geotechnical
Journal, 46(5):513535, 2009.
[13] Knut H Andersen and Rolf Lauritzsen. Bearing capacity for foundations with cyclic
loads. Journal of Geotechnical Engineering, 114(5):540555, 1988.
[14] Knut H Andersen, Arne Kleven, and Dag Heien. Cyclic soil data for design of gravity
structures. Journal of Geotechnical Engineering, 114(5):517539, 1988.
[16] CC Ladd. Stress-deformation and strength characteristics, state of the art report.
Proc. of 9th ISFMFE., 1977, 4:421494, 1977.
[17] S Gourvenec and S Barnett. Undrained failure envelope for skirted foundations under
general loading. Geotechnique, 61(3):263, 2011.
[18] HP Jostad, F Nadim, and KH Andersen. A computational model for fixity of spud
cans on stiff clay. In International Journal of Rock Mechanics and Mining Sciences
and Geomechanics Abstracts, volume 32. Elsevier Science, 1995.
Appendix A
47
48 Soil normalized cyclic shear strength
0 E /sC
cy,f 0 0.48 0.48 0 0.52 0.52
10 u
DSS /sDSS
cy,f 0 1 1 0 1 1
u
C /sC
cy,f 0 1 1 0 1 1
u
E /sC
cy,f 0 0.48 0.48 0 0.52 0.52
100 u
DSS /sDSS
cy,f 0 1 1 0 1 1
u
C /sC
cy,f 0.25 0.95 1.2 0.52 0.76 1.28
u
E /sC
cy,f 0.18 0.38 0.56 0.2 0.46 0.66
1 u
DSS /sDSS
cy,f 0.44 0.87 1.31 0.47 0.94 1.41
u
C /sC
cy,f 0.23 0.93 1.16 0.44 0.6 1.04
u
0.5 E /sC
cy,f 0.17 0.31 0.48 0.15 0.43 0.58
10 u
DSS /sDSS
cy,f 0.4 0.8 1.2 0.38 0.77 1.15
u
C /sC
cy,f 0.22 0.82 1.04 0.32 0.44 0.76
u
E /sC
cy,f 0.16 0.25 0.41 0.1 0.36 0.46
100 u
DSS /sDSS
cy,f 0.35 0.7 1.05 0.27 0.54 0.81
u
49
1 E /sC
cy,f 0.3 0.1 0.4 0.22 0.37 0.59
10 u
DSS /sDSS
cy,f 0.64 0.64 1.28 0.58 0.58 1.16
u
C /sC
cy,f 0.36 0.58 0.94 0.33 0.3 0.63
u
E /sC
cy,f 0.2 0.2 0.4 0.08 0.33 0.41
100 u
DSS /sDSS
cy,f 0.51 0.51 1.02 0.36 0.36 0.72
u
C /sC
cy,f 0.84 0.42 1.26 0.72 0.16 0.88
u
E /sC
cy,f 0.57 0.13 0.7 0.61 0.12 0.73
1 u
DSS /sDSS
cy,f 1.16 0.58 1.74 1 0.5 1.5
u
C /sC
cy,f 0.56 0.42 0.98 0.52 0.15 0.67
u
2 E /sC
cy,f 0.4 0.007 0.407 0.42 0.1 0.52
10 u
DSS /sDSS
cy,f 0.84 0.42 1.26 0.74 0.37 1.11
u
C /sC
cy,f 0.43 0.43 0.86 0.24 0.1 0.34
u
E /sC
cy,f 0.3 0.05 0.35 0.3 0.09 0.39
100 u
DSS /sDSS
cy,f 0.66 0.34 1 0.43 0.23 0.66
u
C /sC
cy,f 0.8 0 0.8 0.67 0 0.67
u
E /sC
cy,f 0.8 0 0.8 0.67 0 0.67
1 u
DSS /sDSS
cy,f 1.28 0 1.28 1.1 0 1.1
u
C /sC
cy,f 0.55 0 0.55 0.47 0 0.47
u
a = 0 E /sC
cy,f 0.55 0 0.55 0.47 0 0.47
10 u
DSS /sDSS
cy,f 0.92 0 0.92 0.78 0 0.78
u
C /sC
cy,f 0.4 0 0.4 0.28 0 0.28
u
E /sC
cy,f 0.4 0 0.4 0.28 0 0.28
100 u
DSS /sDSS
cy,f 0.69 0 0.69 0.45 0 0.45
u
50 Soil normalized cyclic shear strength
Appendix B
51
52 Soil stress strain curves
1.2 1.2
1 1
0.8 0.8
0.6 0.6
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.4 0.4
0.2 0.2
0 0
-20 -10 -0.2 0 10 20 -30 -20 -10 -0.2 0 10 20 30 40
-0.4 -0.4
-0.6 -0.6
-0.8 -0.8
(%) (%)
(a) cy = 0 (b) cy /a = 1
1.2 1
1 0.8
0.8
0.6
0.6
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.4 0.4
0.2 0.2
0 0
-20 -10 -0.2 0 10 20 -20 -10 0 10 20
-0.2
-0.4
-0.6 -0.4
-0.8 -0.6
(%) (%)
(c) cy /a = 2 (d) a = 0
1.2 2
1
1.5
0.8
0.6 1
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.4
0.5
0.2
0 0
-20 -10 0 10 20 -20 -10 0 10 20
-0.2
-0.5
-0.4
-0.6 -1
(%) (%)
(a) cy = 0 (b) cy /a = 1
2 1.6
1.4
1.5
1.2
1
1
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.8
0.5 0.6
0.4
0
0.2
-20 -10 0 10 20
0
-0.5
-20 -10 -0.2 0 10 20
-1 -0.4
(%) (%)
(c) cy /a = 2 (d) a = 0
1.2 1.4
1 1.2
0.8 1
0.8
0.6
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.6
0.4
0.4
0.2
0.2
0
0
-20 -10 0 10 20
-0.2 -20 -10 -0.2 0 10 20
-0.4 -0.4
-0.6 -0.6
(%) (%)
(a) cy = 0 (b) cy /a = 1
1.4 1
1.2
0.8
1
0.6
0.8
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.6 0.4
0.4 0.2
0.2
0
0 -20 -10 0 10 20
-20 -10 0 10 20 -0.2
-0.2
-0.4 -0.4
(%) (%)
(c) cy /a = 2 (d) a = 0
1.2 1.2
1 1
0.8
0.8
0.6
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.4 0.6
0.2 0.4
0 0.2
-20 -10 -0.2 0 10 20
0
-0.4 -20 -10 0 10 20
-0.6 -0.2
-0.8 -0.4
(%) (%)
(a) cy = 0 (b) cy /a = 1
1.2 0.9
0.8
1
0.7
0.8 0.6
cy/suC or cy/suDSS
cy/suC or cy/suDSS
0.5
0.6 0.4
0.4 0.3
0.2
0.2 0.1
0
0
-20 -10 -0.1 0 10 20
-20 -10 0 10 20
-0.2 -0.2
(%) (%)
(c) cy /a = 2 (d) a = 0
MATLAB code
57
58 MATLAB code
%% The MATLAB code in this section was made for linear interpolation bearing
capacity analysis results for general loading condition to plot the cross-
section failure envelopes in 3D failure surface.
clear all
clc
close all
%% plot three loading combined situation
n M=0.2:0.1:0.90;
cols n M=size(n M,2);
% caculate the H loads for M/M max=0.4:0.1:0.7 without vertical loading
N 100data=xlsread('C:\Users\yo\Dropbox\NTNU\Master Thesis\HVMCap for matlab
programming\Read data.xlsx','sheet5');
% HMdata=N 100data(:,(8*(14-1)+1):(8*(18-1)+8));
HMdata=N 100data(:,(8*(14-1)+1):(8*(25-1)+8));
MVdata=N 100data(:,(8*(4-1)+1):(8*(8-1)+8));
HVdata=N 100data(:,(8*(9-1)+1):(8*(13-1)+8));
cols HMdata=40;% coloumns of the HMdata sets
for i=1:5 % caculate the maximum HVM loads for different combinations for HM&
MV plane
C.2 MATLAB code for plotting 3D failure surface 61
for i=1:7
for j=1:5
H N1 re(j,i)=H N1(5*(i-1)+j);
V N1 re(j,i)=V N1(5*(i-1)+j);
M N1 re(j,i)=M N1(5*(i-1)+j);
end
end
zero=zeros(1,cols n M);% for compensate the 0 values for V load in HM plane
and H load in MV Plane
one=ones(5,1);% for compensate the 1 values for pure moment loading
zerore=zeros(5,1);
rows M N1=size(M N1 re ,1);
for i=1:cols n M
for j=1:rows M N1
H N1VsM(j,i)=interp1(M N1 re(j,:),H N1 re(j,:),n M(i));
V N1VsM(j,i)=interp1(M N1 re(j,:),V N1 re(j,:),n M(i));
end
end
one=ones(7,1);% for compensate the 1 values for pure moment loading
zerore=zeros(7,1);
H loadV re=[1 H loadV 0];% h load only in HV plane M/Mmax=0
V loadH re=[0 V loadH 1];% v load only in HV plane M/Mmax=0
H N1VsM re=[H loadVsM;H N1VsM;zero];
V N1VsM re=[zero;V N1VsM;V loadVsM];
H N1VsM ne=[(H loadV re)' H N1VsM re zerore];
V N1VsM ne=[(V loadH re)' V N1VsM re zerore];
rows M N1=size(M N1 re ,1);
n M zero=zeros(1,1);
n M one=ones(1,1);
n M ne=[n M zero n M n M one];
n M re=[n M ne;n M ne;n M ne;n M ne;n M ne;n M ne;n M ne];
% plot HV load combination for constant moment loading
figure(1)
hold on
plot3(H N1VsM ne(:,1),V N1VsM ne(:,1),n M re(:,1),'--rs');
plot3(H N1VsM ne(:,cols n M+2),V N1VsM ne(:,cols n M+2),n M re(:,cols n M+2),'
C.3 MATLAB code for plotting 3D displacement contour diagrams 63
--gs');
for v=2:(cols n M+1)% suppose 1:cols n M+2
plot3(H N1VsM ne(:,v),V N1VsM ne(:,v),n M re(:,v),'--bs');
xlabel('H/H u l t');
ylabel('V/V u l t');
zlabel('M/M u l t');
legend('M/M u l t=0','M/M u l t=1','M/M u l t=0.3 - 0.7','Location','NorthWest
');
legend('boxoff');
grid on
diagrams
%% The following MATLAB script was made for plotting displacement for combined
horizontal, vertical and overturning moment load
clear all
clc
close all
%% load 3D combined load situations results into workspace
% define the interval of rotational rad
n theta=0.001:0.001:0.14;% rotational rad
cols n theta=size(n theta,2);% cols of n M
n comb=35;%number of load combinations
n M=0.7;%:0.1:1;%M/M max=0.3,0.4,0.5
cols n M=size(n M,2);%cols of n M
% read data
N 100data=xlsread('C:\Users\yo\Dropbox\NTNU\Master Thesis\HVMCap for matlab
programming\Read data.xlsx','sheet10');
depth=xlsread('C:\Users\yo\Dropbox\NTNU\Master Thesis\HVMCap for matlab
programming\Load combinations allPlane Geo.1.xlsx','sheet1','E2');% bucket
foundation height
% read the M ult,V ult,H ult
M ult=max(N 100data(:,8*(1-1)+6));
H ult=max(N 100data(:,8*(2-1)+7));
64 MATLAB code
V ult=max(abs(N 100data(:,8*(3-1)+8)));
% interpolate each data set for different rotational rad
interp data=zeros(35* cols n theta,7);% workspace for putting interpolate
results
for i=4:38
x Tdisp=N 100data(:,8*(i-1)+2);% top hori. disp.
x Tdisp(x Tdisp==0)=[];% delete zero
x Bdisp=N 100data(:,8*(i-1)+4);% bottom hori. disp.
x Bdisp(x Bdisp==0)=[];% delete zero
z Tdisp=N 100data(:,8*(i-1)+3);% top vert. disp.
z Tdisp(z Tdisp==0)=[];% delete zero
z Bdisp=N 100data(:,8*(i-1)+5);% bottom vert. disp.
z Bdisp(z Bdisp==0)=[];% delete zero
M=N 100data(:,8*(i-1)+6);% moment load
M(M==0)=[];% delete zero
M max=max(abs(M));% read the maximum moment for each load combination
H=N 100data(:,8*(i-1)+7);% hori. load
H(H==0)=[];% delete zero
H max=max(abs(H));% read the maximum moment for each load combination
V=N 100data(:,8*(i-1)+8);% Vertical load
V(V==0)=[];% delete zero
V max=max(abs(V));% read the maximum moment for each load combination
%calculate theta;
H disp=0.67* x Bdisp+0.33* x Tdisp;%calculate horizontal disp
theta=atan((x Tdisp-H disp)./(2/3*depth)); %calculate rotation rad
[thetaUni,ia]=unique(theta,'stable');% monotonic increase
MUni=M(ia);% pick out moment according to monotonic increase theta
HUni=H(ia);% pick out horizontal load according to monotonic increase theta
VUni=V(ia);% pick out vertical load according to monotonic increase theta
interp data(1:cols n theta,7*(i-4)+1)=n theta;% read correspoding maximum
vertical load for each load combinations
interp data(1:cols n theta,7*(i-4)+2)=interp1(thetaUni,MUni,n theta);%
interpolate theta to get moment load
interp data(1:cols n theta,7*(i-4)+3)=interp1(thetaUni,HUni,n theta);%
interpolate theta to get moment load
interp data(1:cols n theta,7*(i-4)+4)=interp1(thetaUni,VUni,n theta);%
C.3 MATLAB code for plotting 3D displacement contour diagrams 65
H fi(7*(k-1)+j,i)=H new(5*(j-1)+k,i);
V fi(7*(k-1)+j,i)=V new(5*(j-1)+k,i);
M max fi(7*(k-1)+j,i)=M max new(5*(j-1)+k,i);
H max fi(7*(k-1)+j,i)=H max new(5*(j-1)+k,i);
V max fi(7*(k-1)+j,i)=V max new(5*(j-1)+k,i);
end
end
% interpolate the M/M max to get different H,V for the same theta
for m=1:5
H fi new(1:7,1)=H fi(((7*(m-1)+1):(7*(m-1)+7)),i);
V fi new(1:7,1)=V fi(((7*(m-1)+1):(7*(m-1)+7)),i);
M max fi new(1:7,1)=M max fi(((7*(m-1)+1):(7*(m-1)+7)),i);
H max fi new(1:7,1)=H max fi(((7*(m-1)+1):(7*(m-1)+7)),i);
V max fi new(1:7,1)=V max fi(((7*(m-1)+1):(7*(m-1)+7)),i);
[H fi new,ia]=unique(H fi new,'stable');
V fi new=V fi new(ia);
M max fi new=M max fi new(ia);
H fi new re(m,i)=interp1(M max fi new ,H fi new,n M);
V fi new re(m,i)=interp1(M max fi new ,V fi new,n M);
H max fi new re(m,i)=interp1(M max fi new ,H max fi new ,n M);
V max fi new re(m,i)=interp1(M max fi new ,V max fi new ,n M);
end
end
%% select H & V for different M/M max and theta in HM and VM plan,seperately
%read data in HM plane and VM plane
two plane data=xlsread('C:\Users\yo\Dropbox\NTNU\Master Thesis\HVMCap for
matlab programming\Read data.xlsx','sheet9');
HMdata=two plane data(:,(8*(14-1)+1):(8*(18-1)+8));
MVdata=two plane data(:,(8*(4-1)+1):(8*(8-1)+8));
cols HMdata=40;% coloumns of the load plane data sets
%% interpolate theta to get different H,M, H max & M max in HM plane
for i=1:cols HMdata/8% 2D load plane data sets
x Tdisp=HMdata(:,8*(i-1)+2);% top hori. disp.
x Tdisp(x Tdisp==0)=[];% delete zero
x Bdisp=HMdata(:,8*(i-1)+4);% bottom hori. disp.
C.3 MATLAB code for plotting 3D displacement contour diagrams 67
end
end
% divide different theta into different conlumns
for i=1:cols n theta
theta newHM(:,i)=theta reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+
n 2loadcomb));
H newHM(:,i)=H reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+n 2loadcomb))
;
M newHM(:,i)=M reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+n 2loadcomb))
;
M max newHM(:,i)=M max reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+
n 2loadcomb));
H max newHM(:,i)=H max reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+
n 2loadcomb));
end
% interpolate M/M max in HM plane
for i=1:cols n theta;
[H fi,ia]=unique(H newHM(:,i),'stable');% delete
M max newHM re=M max newHM(:,i);
M max HM=M max newHM re(ia);
H HM(cols n M,i)=interp1(M max HM,H fi,n M);
end
%% interpolate theta to get different V,M, V max & M max in VM planeHdisp no
this part
for i=1:cols HMdata/8% 2D load plane data sets
x Tdisp=MVdata(:,8*(i-1)+2);% top hori. disp.
x Tdisp(x Tdisp==0)=[];% delete zero
x Bdisp=MVdata(:,8*(i-1)+4);% bottom hori. disp.
x Bdisp(x Bdisp==0)=[];% delete zero
z Tdisp=MVdata(:,8*(i-1)+3);% top vert. disp.
z Tdisp(z Tdisp==0)=[];% delete zero
z Bdisp=MVdata(:,8*(i-1)+5);% bottom vert. disp.
z Bdisp(z Bdisp==0)=[];% delete zero
M=MVdata(:,8*(i-1)+6);% moment load
M(M==0)=[];% delete zero
M max=max(abs(M));% read the maximum moment for each load combination
C.3 MATLAB code for plotting 3D displacement contour diagrams 69
;
M max newHM(:,i)=M max reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+
n 2loadcomb));
V max newHM(:,i)=V max reHM((n 2loadcomb*(i-1)+1):(n 2loadcomb*(i-1)+
n 2loadcomb));
end
% interpolate M/M max in MV plane
for i=1:cols n theta;
[V fi,ia]=unique(V newHM(:,i),'stable');% delete
M max newHM re=M max newHM(:,i);
M max HM=M max newHM re(ia);
V HM(cols n M,i)=interp1(M max HM,V fi,n M);
end
%% plot the displacement contours for a specfic M/M max
zero=zeros(1,cols n theta);% compensate the zeros for MV for H=0 and HM for V
=0
H plot=[H HM;H fi new re;zero];% integarate 2load plane data into 3 load plane
data
V plot=[zero;V fi new re;V HM];
figure(1)
hold on
for i=1:cols n theta
x=H plot(:,i);
y=V plot(:,i);
plot(x,y,'-k*');
axis([0 1 0 1]);
end
title('M/M u l t=0.7');
xlabel('H/H u l t');
ylabel('V/V u l t');
legend('\theta=0.001 - 0.14 rad','Location','Best');%,'k M');
legend('boxoff');
axis([0 1 0 1]);
grid on
Appendix D
H=0 M=0
1 1
= 0.0010.14 rad u = 0.011 m
0.9 w = 0.011 m 0.9 w = 0.011 m
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
71
72 Displacement contour diagrams
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
u = 0.011 m
w = 0.011 m
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
0.9 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
= 0.0010.14 rad
0.9 w = 0.011 m 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
0.1 u = 0.011 m
= 0.0010.14 rad
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
H=0 M=0
1 1
0.9 0.9
0.8 0.8
0.7 0.7
0.6 0.6
M/Mult
H/Hult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
= 0.0010.14 rad u = 0.011 m
w = 0.011 m w = 0.011 m
0.1 0.1
Failure envelope Failure envelope
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
V/Vult V/Vult
V=0
1
0.9
0.8
0.7
0.6
M/Mult
0.5
0.4
0.3
0.2
u = 0.011 m
0.1 = 0.0010.14 rad
Failure envelope
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
M/Mult=0.6 M/Mult=0.6
1 1
0.9 0.9
=0.001 0.14 rad u = 0.01 1 m
0.8 0.8
0.7 0.7
0.6 0.6
V/Vult
V/Vult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult H/Hult
M/Mult=0.6
1
0.9
w = 0.01 1 m
0.8
0.7
0.6
V/Vult
0.5
0.4
0.3
0.2
0.1
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult
M/Mult=0.7 M/Mult=0.7
1 1
0.9 0.9
=0.001 0.14 rad u = 0.01 1 m
0.8 0.8
0.7 0.7
0.6 0.6
V/Vult
V/Vult
0.5 0.5
0.4 0.4
0.3 0.3
0.2 0.2
0.1 0.1
0 0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult H/Hult
M/Mult=0.7
1
0.9
w = 0.01 1 m
0.8
0.7
0.6
V/Vult
0.5
0.4
0.3
0.2
0.1
0
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
H/Hult