Energy Calculation
Introduction
Energy changes determined in 3DEC are performed for the intact rock, the joints, and for the work done on boundaries. The energy terms calculated here use the same general nomenclature as those used by Salamon (1984).
Since 3DEC uses an incremental solution procedure, the equations of motion are solved at each mass point in the body at every timestep. The incremental change in energy components is determined at each timestep as the system attempts to come to equilibrium. 3DEC also keeps a running sum of each component.
Energy Balance
The total energy balance can be expressed in terms of the released energy (\(W_r\)), which is the difference between the work done at the boundary of the model and the total stored and dissipated strain energies:
where:
\(Wr\) = released energy;
\(W\) = total boundary loading work supplied to the system;
\(U_c\) = total stored strain energy in material;
\(U_b\) = total change in potential energy of the system;
\(W_j\) = total dissipated energy in joint shear; and
\(W_p\) = total dissipated work in plastic deformation of intact rock.
A second calculation of released energy can be made based on the kinetic energy, mass damping work, the work performed at viscous boundaries, and the strain energy in excavated material:
where:
\(U_k\) = current value of kinetic energy in the system;
\(W_k\) = total work dissipated by mass damping;
\(W_v\) = work done by viscous (non-reflecting) boundaries; and
\(U_m\) = total strain energy in excavated material.
This second form of the released energy is particularly useful for dynamic problems because the released kinetic energy is easily calculated.
The definitions of these individual terms, and the method of their calculation, are given in the following section.
Calculation of Individual Energy Components
Total Boundary Loading Work (\(W\))
The work done on the external boundaries of the 3DEC block structure is calculated from the boundary gridpoint forces and displacements. For dynamic problems, the stress boundary may be replaced by non-reflecting (viscous) boundaries after initial stress equilibrium is achieved. The work done at a gridpoint is calculated:
where:
\(W_{gi}\) = work done at gridpoint \(i\) on the outer boundary;
\(F_{xi}\) = \(x\)-oriented force at the gridpoint;
\(u_{xi}\) = \(x\)-oriented displacement at the gridpoint;
\(F_{yi}\) = \(y\)-oriented force at the gridpoint; and
\(u_{yi}\) = \(y\)-oriented displacement at the gridpoint.
\(F_{zi}\) = \(z\)-oriented force at the gridpoint; and
\(u_{zi}\) = \(z\)-oriented displacement at the gridpoint.
The total boundary work (\(W_l\)) done during a timestep (\(l\)) is the sum of the work done for all gridpoints on the boundary:
where \(ngp\) = the number of gridpoints on the boundary.
The boundary work is summed for all timesteps:
where \(nt\) = the number of timesteps. \(W\) will approach a constant value as the system approaches equilibrium.
The rate of boundary loading work, \(\Delta W\), is simply the boundary work done per timestep:
where \(\Delta t\) is the timestep.
Histories of the total work, \(W\), and incremental work, \(\Delta W\), are kept during a simulation.
Potential Energy (\(U_b\))
The change in gravitational potential energy is calculated from the gridpoint gravitational forces and the displacements of the gridpoints. The total potential energy is summed for incremental displacements of all gridpoints:
where:
\(U_{gi}\) = the potential energy of gridpoint \(i\);
\(m_i\) = the mass of gridpoint \(i\);
\(u_{xi}\), \(u_{yi}\), \(u_{zi}\) = displacement components of gridpoint \(i\); and
\(g_x\), \(g_y\), \(g_z\) = accelerations in the \(x\)-, \(y\)-, and \(z\)-directions (usually gravity).
The total potential energy is found by summing the energy for all gridpoints, \(i\), at a given timestep, \(k\):
The total \(U_b\) is kept for all timesteps:
where \(nt\) is the number of timesteps; \(U_b\) will approach a constant value as the model approaches equilibrium.
Kinetic Energy (\(U_k\))
The kinetic energy is determined for each gridpoint at each timestep, and is summed for all gridpoints at that timestep. A running total of the kinetic energy is not kept; so, as the system approaches equilibrium, the kinetic energy will approach zero. The kinetic energy is given by
where:
\(U_k\) = kinetic energy of all gridpoints in a given timestep;
\(m_i\) = mass of gridpoint \(i\); and
\(\dot u_i\) = velocity at gridpoint \(i\).
The incremental kinetic energy is calculated so the user can examine the rate of change between timesteps:
where:
\(U_{kn}\) = kinetic energy at timestep \(n\);
\(U_{kn−1}\) = kinetic energy at timestep \(n-1\); and
\(\Delta t\) = timestep.
Kinetic energy is also calculated for rigid blocks in a similar fashion, based on the block mass and velocity.
Damped Energy (\(W_k\))
The mass-damping work is the summation of all energy absorbed by either local damping or global damping, and is intended for use primarily with static analysis. For dynamic analyses, the work done on viscous boundaries will be much larger than the damped energy, and will largely control the calculated value of the total released energy.
The damped energy can most easily be seen by examining a simplified version of the equation of motion,
where:
\(\dot u\) = velocity of a gridpoint of mass, \(m\)
\(\sum F\) = the force sum at the gridpoint; and
\(\alpha\) = the damping coefficient, which is \(2 \pi f \gamma\) , where \(\gamma\) = fraction of critical damping; and \(f\) = natural frequency of the system (cps); the circular frequency in radians per second is \(\omega = 2 \pi f\).
The damping force is given by
and the rate of damped energy change at a gridpoint is
The damped energy over a timestep at a gridpoint, \(j\), is
where \(W_{dj}\) is the energy damped at gridpoint \(j\), and \(U_k\) is the kinetic energy of the gridpoint.
Therefore, the damped and kinetic energy components are related by the damping coefficient. The total mass damping work is the sum of all gridpoints and timesteps:
where:
\(W_d\) = total energy damped;
\(nt\) = number of timesteps; and
\(ngp\) = number of gridpoints.
The incremental damped energy is also calculated as
Strain Energy Stored in the Rock Mass (\(U_c\))
The total strain energy stored in the rock mass is composed of two parts: the energy stored in the blocks, and that stored in the joints. Each is calculated, and the total stored energy, \(U_c\), is determined as the sum of these two components.
Block-Stored Strain Energy (\(U_{cb}\))
The strain energy in the blocks is determined for all finite difference zones during each timestep. The total stored strain energy is calculated by summing the values for all blocks. The incremental strain energy in each zone for a timestep is given by
where:
\(\sigma_{ij}\) = current zone stresses;
\(\sigma_{ij}^\prime\) = zone stresses from the previous timestep;
\(e_{ij}\) = incremental strains over the current timestep; and
\(V\) = volume of zone.
The strain energy in a block is the sum of all zones within the block,
The total strain energy in a given timestep is the sum for all blocks. A running total is kept for all timesteps, so that the value of \(U_{c_b}\) will approach a constant value with time as the system approaches equilibrium.
The incremental block strain energy is also calculated by
An exception to the above calculation method is used in the case of initial equilibrium of the model prior to any excavation. Normally, the in-situ stresses are set up in the model by “freezing” the stresses in each zone using the block insitu stress command. By also applying boundary stresses (which are in equilibrium with the internal stresses) at the same time, the model will be at initial stress equilibrium and will not require any timestepping. However, this method of consolidating the body under initial stresses produces no strain in the body, and no apparent strain energy. Therefore, the alternate form of the strain-energy density equation is used to define the initial strain energy in the system:
Joint Strain Energy (\(U_{cj}\))
The strain energy stored in the joints is separated into four component parts for elastic strain in shear (\(U_{js}\)), compression (\(U_{jc}\)), tension (\(U_{jt}\)) and energy dissipated in slip (\(U_{jf}\)). The calculations of these components are performed differently for the Mohr-Coulomb-based constitutive models and the Continuously Yielding model. In the Mohr-Coulomb-based models, the joint normal and shear stiffnesses are linear, whereas the continuously yielding model allows nonlinear stiffness. The energy is determined for each contact along all joints in the model, although, at present, it is not possible to separate energy by joint.
Coulomb Slip Model, Linear Stiffness
If \(f_n < 0\),
If \(f_n \ge 0\),
If \(f_s < f_{s \max}\),
If \(f_s \ge f_{s \max}\),
where:
\(f_n\), \(f_s\) = current normal shear force at a contact, compression positive;
\(f_n^{\prime}\), \(f_s^{\prime}\) = previous normal shear forces at a contact;
\(u_n\), \(u_s\) = incremental normal and shear displacements at the contact over the current timestep; and
\(f_{s \max}\) = shear stress at which the Coulomb slip condition is met (\(f_{s \max} < f_n \tan \phi + C\)).
For the case in which \(f_s \leq f_{s \max}\) and slip occurs, energy is dissipated in the Coulomb model by friction (see below).
Continuously Yielding Model – Nonlinear Stiffness Allowed
where:
\(f_n\), \(f_s\) = normal shear forces at a contact;
\(f_n^{\prime}\), \(f_s^{\prime}\) = previous contact stresses;
\(u_n\), \(u_s\) = incremental, normal shear displacements over the current timestep; and
\(F\) = the yielding factor for the continuously yielding model as described in i Continuously Yielding Joint Model In 3DEC.
The total energy absorbed and dissipated by the joints, \(U_{cj}\), is given by the sum of all components:
Strain Energy Content of Excavated Material (\(U_m\))
When rock is excavated, the strain energy that was stored in the excavated volume is released. 3DEC allows excavation of the rock blocks that form the opening through use of the block delete command or assignment of the null constitutive model (via block excavate). The null model does not delete the blocks, but forces in null blocks do not pass to gridpoints of adjoining blocks. The null zones can collapse due to deformation of the opening, and can later be changed to a backfill material. When a block is deleted or given a null constitutive model, the energy sums are updated.
The total strain energy in the excavated material consists of the strain energy in the blocks and the joints. The strain energy in the blocks is calculated in the same manner as the strain energy described previously. Using principal stresses this can be calculated by:
where:
\(U_{mb}\) = block strain energy;
\(V\) = volume of zone;
\(E\) = Young’s modulus of the rock mass;
\(\nu\) = Poisson’s ratio;
\(\sigma_{1}\), \(\sigma_{2}\), \(\sigma_{3}\) = principal stresses in zone centroid;
\(nb\) = number of blocks in excavated material; and
\(nz\) = number of zones in the block.
The total strain energy in the joints bounding the excavated block(s) (\(U_{cj}\)) is given as follows.
Coulomb Joints (Constant Stiffness)
where:
\(U_{cj}\) = strain energy stored in the joints;
\(f_n\), \(f_s\) = normal shear force in joints;
\(k_n\), \(k_s\) = normal shear stiffness of joints; and
\(nc\) = number of contacts.
Continuously Yielding Joints (Nonlinear Normal Stiffness)
where:
\(e_n\) = normal stiffness exponent;
\(u_n\) = normal displacement (closure) of joint surfaces;
\(a_n\) = initial normal stiffness of joint;
\(f_s\) = shear force at contact;
\(k_s\) = shear stiffness of contact; and
\(nc\) = number of contacts.
When a block is excavated, its energy is removed from the total strain energy, \(U_c\), and added to the total for the excavated material, \(U_{mb}\). The initial or in-situ strain energy state for the rock mass is determined by using the standard strain-energy density function, where the principal stresses are equal to the in-situ stresses. This is added to the boundary loading work (W) for the initial equilibrium, pre-mining condition. The final values for stored strain energies are determined:
where a “ \(^\prime\) ” refers to the values from previous excavation steps.
Friction Work Done on Joints (\(W_j\))
Energy is dissipated through frictional heating of joints. This work done is exchanged from the elastic strain energy, and is irrecoverable. 3DEC keeps track of the frictional energy separately from the elastic (stored) joint energy terms (\(U_{jt}\), \(U_{jc}\) and \(U_{js}\)). The friction loss is calculated for linear and nonlinear normal stiffness as follows.
Coulomb Joints, Linear Normal Stiffness
If \(f_s \ge f_{s \max}\),
where:
\(U_{jf}\) = frictional energy at the contact during a timestep;
\(f_s\) = current shear force at a contact;
\(f_s^{\prime}\) = previous shear force at a contact;
\(u_s\) = increment in shear displacement; and
\(nc\) = number of contacts.
Continuously Yielding Joint (Nonlinear Normal Stiffness)
where:
\(F\) = the yield factor for the continuously yielding joint; and
\(nc\) = number of contacts.
The total dissipated energy is kept by summing over all the timesteps during an excavation step,
where:
\(W_j\) = total dissipated friction energy; and
\(nt\) = number of timesteps.
Viscous Boundary Work (\(W_v\))
Viscous boundaries are used to dampen reflections of incident stress waves. The energy damped from these stress waves is calculated from the boundary forces and deflections at the boundary gridpoints:
where:
\(W_{gj}\) = boundary work at a gridpoint, \(j\);
\(f_x\), \(f_y\), \(f_z\) = boundary forces; and
\(u_x\), \(u_y\), \(u_z\) = boundary displacements.
The viscous energy for a timestep is given by
where \(nbp\) = number of boundary gridpoints.
The total viscous work is summed for all timesteps. The incremental viscous boundary work is calculated by
where \(\Delta t\) = the timestep.
Energy Dissipation in Blocks through Plastic Work (\(W_p\))
Several plasticity models that can describe the deformability of the blocks are available in 3DEC. Energy is dissipated through plastic work as the zones undergo irreversible deformation. The strain in any zone can be divided into an elastic and a plastic part. The elastic strain can be determined, followed by the elastic strain energy as determined previously. The plastic work is found by taking the difference between the total strain energy and the elastic energy component.
The elastic strain energy is given by
The change in elastic strain energy between excavation steps is given by
where \(W_e\) is the current elastic strain energy, and \(W_e^{\prime}\) is the previous elastic strain energy.
The total energy change can be found from the total strain and stress
The total plastic work dissipated during an excavation step is the difference between the total and elastic energy change at any timestep,
and the total dissipated is simply the sum of \(\Delta W_p\) for all blocks at each timestep.
Energy Dissipated in Backfill Compression
3DEC keeps track of those blocks that have been excavated using the null model method and replaced by backfill. In this case, elastic strain energy may be stored in the fill, and energy may be dissipated via plastic deformation. This would normally be the case for sandfills. Energy values are calculated as either stored or dissipated for each zone (as previously described), and are added to the plastic work term.
Volume of Excavated Material (\(V_m\))
When a block is deleted or assigned a null constitutive model type, the volume is added to the value \(V_m\). The volume is calculated by using the mass and density of the zones that compose it:
\(V_z\) = mass/density
where \(V_z\) = volume of a zone.
The volume of the deleted blocks is then equal to
for all zones in the excavated blocks.
Method of Operation in 3DEC
The energy calculations in 3DEC are initiated using the model config energy command. From this point on, the energies described in the previous sections are calculated in an incremental fashion at each timestep from the stress, force, displacement and strain changes. All energy values are summed from this point, with the exception of the kinetic energy, \(U_k\), which is kept as an incremental value. Therefore, the magnitude of the energies upon printout will be the sum for the problem since timestepping began, and will include that computed for all excavation steps.
Energy quantities can be monitored using model history energy commands. At any time, the energy values can be printed to the screen with the block list energy command.
The following section presents an example.
Example: Excavation of a Circular Hole
The model shows the evolution of energy components after rapid excavation of a circular tunnel. The radius of the excavation is 1 m. The outer radius of the 3DEC model is 10 m. A hydrostatic compressive stress of 100 MPa exists prior to excavation. The elastic material has a shear modulus of 29.17 GPa and Poisson’s ratio of 0.2.
Figure 1 shows the 3DEC block geometry with zones. The 3DEC data file to simulate the excavation and monitor energy components is listed in energy.dat
Figure 1: 3DEC model geometry
energy.dat
model new
model random 10000
model configure energy
model large-strain off
block create drum center-1 0 0 0 center-2 0 1 0 ...
radius-1 10 radius-2 10 edges 24
; make table to define tunnel coordinates
; input table name, radius and number of edges
fish define make_table(table_name, radius, num_edges)
tab = table.create(table_name)
loop i (1,num_edges)
local theta = 2.0*math.pi*(i-1)/num_edges
table.x(tab,i) = radius*math.cos(theta)
table.y(tab,i) = radius*math.sin(theta)
end_loop
end
[make_table('tunnel',1,24)]
block cut tunnel radial table-1 'tunnel' ...
table-2 'tunnel' axis 0,0,0 0,1,0 group 'tunnel'
block zone size center 0 0 0 edgelength-1 0.2 distance-1 2 ...
edgelength-2 1.5 distance-2 10
block zone generate-new minimum-edge 0.2 maximum-edge 1.5
;material properties (MPa) - Poisson's ratio = 0.2
block zone cmodel assign elastic
block zone property density .002 bulk 38.9e3 shear 29.17e3
;set stresses
block insitu stress -100 -40 -100 0 0 0
; fix boundaries
block gridpoint apply velocity 0 0 0 ...
range cylinder end-1 0 0 0 end-2 0 1 0 radius 9.5 10.5
block gridpoint apply velocity-y 0 range position-y 0
block gridpoint apply velocity-y 0 range position-y 1
block history displacement-x position 1 0.5 0
block history displacement-x position 2 0.5 0
block history stress-xx position 1 0.5 0
block history stress-xx position 2 0.5 0
block history stress-xx position 3 0.5 0
block history stress-zz position 1 0.5 0
block history stress-zz position 2 0.5 0
block history stress-zz position 3 0.5 0
model solve convergence 1
model save 'initial'
============================================
; excavation
model history energy block-kinetic
model history energy block-strain
model history energy block-strain-excavated
model history energy block-volume-excavated
model history energy damp-mass
model history energy boundary-load
model history energy total-released
block delete range group 'tunnel'
model solve convergence 1
program log-file 'energies.log'
program log on
block list energy
program log off
model save 'exc'
Figure 2 shows a history of the strain, kinetic and damping energies. In this figure, the kinetic energy term is not summed over time, but decays to zero as the model comes to equilibrium. At the same time, damped energy, which is summed, approaches a constant value. Similarly, the strain energy approaches a constant negative value — strain energy is released after the tunnel is excavated and the rock moves inwards.
Figure 2: History of the damped, strain and incremental kinetic energies.
The table below shows the result of the block list energy command.
**********************************************
* Logging started at Tue Dec 17 00:08:46 2024
* By 3dec Version 9.3 Release 24
*
*
* Job Title:
**********************************************
Totals for energy stored and dissipated in system.
--------------------------------------------------
Current kinetic energy (Uk) = 4.336171e-05
Total block strain energy (Ucb) = 3.110524e+01
Total fill strain energy (Ucf) = 0.000000e+00
Total joint strain energy (Ucj) = 0.000000e+00
Total strain energy (Uc=Ucb+Ucf+Ucj) = 3.110524e+01
Total block energy excavated (Umb) = 3.193760e-01
Total joint energy excavated (Umj) = 0.000000e+00
Total strain energy excavated (Um=Umb+Umj) = 3.193760e-01
Total block volume excavated (Vm) = 3.105829e+00
Total change in potential energy (Ub) = 0.000000e+00
Total mass damping work (Wk) = 0.000000e+00
Total viscous boundary work (Wv) = 0.000000e+00
Total friction work (Wj) = 0.000000e+00
Total plastic strain work (Wp) = 0.000000e+00
Total boundary loading work (W) = 3.193760e+01
Total energy stored (Wr=W-Uc-Ub-Wj-Wp) = 8.323569e-01
Total energy released (Wr=Uk+Wk+Wv+Um) = 3.194193e-01
Breakdown of energy stored in joints (Ucj)
------------------------------------------
Total energy stored in tension (Ujt) = 0.000000e+00
Total energy stored in compression (Ujc) = 0.000000e+00
Total energy stored in shear (Ujs) = 0.000000e+00
**********************************************
* Logging ended at Tue Dec 17 00:08:46 2024
**********************************************
**********************************************
* Logging started at Tue Dec 17 22:28:25 2024
* By 3dec Version 9.3 Release 24
*
*
* Job Title:
**********************************************
**********************************************
* Logging ended at Tue Dec 17 22:28:25 2024
**********************************************
**********************************************
* Logging started at Tue Dec 17 23:30:11 2024
* By 3dec Version 9.3 Release 24
*
*
* Job Title:
**********************************************
Totals for energy stored and dissipated in system.
--------------------------------------------------
Current kinetic energy (Uk) = 4.350585e-05
Total block strain energy (Ucb) = 3.110524e+01
Total fill strain energy (Ucf) = 0.000000e+00
Total joint strain energy (Ucj) = 0.000000e+00
Total strain energy (Uc=Ucb+Ucf+Ucj) = 3.110524e+01
Total block energy excavated (Umb) = 3.193760e-01
Total joint energy excavated (Umj) = 0.000000e+00
Total strain energy excavated (Um=Umb+Umj) = 3.193760e-01
Total block volume excavated (Vm) = 3.105829e+00
Total change in potential energy (Ub) = 0.000000e+00
Total mass damping work (Wk) = 0.000000e+00
Total viscous boundary work (Wv) = 0.000000e+00
Total friction work (Wj) = 0.000000e+00
Total plastic strain work (Wp) = 0.000000e+00
Total boundary loading work (W) = 3.193760e+01
Total energy stored (Wr=W-Uc-Ub-Wj-Wp) = 8.323569e-01
Total energy released (Wr=Uk+Wk+Wv+Um) = 3.194195e-01
Breakdown of energy stored in joints (Ucj)
------------------------------------------
Total energy stored in tension (Ujt) = 0.000000e+00
Total energy stored in compression (Ujc) = 0.000000e+00
Total energy stored in shear (Ujs) = 0.000000e+00
**********************************************
* Logging ended at Tue Dec 17 23:30:11 2024
**********************************************
**********************************************
* Logging started at Wed Dec 18 21:25:48 2024
* By 3dec Version 9.3 Release 24
*
*
* Job Title:
**********************************************
**********************************************
* Logging ended at Wed Dec 18 21:25:48 2024
**********************************************
**********************************************
* Logging started at Wed Dec 18 22:27:32 2024
* By 3dec Version 9.3 Release 24
*
*
* Job Title:
**********************************************
Totals for energy stored and dissipated in system.
--------------------------------------------------
Current kinetic energy (Uk) = 4.346831e-05
Total block strain energy (Ucb) = 3.110524e+01
Total fill strain energy (Ucf) = 0.000000e+00
Total joint strain energy (Ucj) = 0.000000e+00
Total strain energy (Uc=Ucb+Ucf+Ucj) = 3.110524e+01
Total block energy excavated (Umb) = 3.193760e-01
Total joint energy excavated (Umj) = 0.000000e+00
Total strain energy excavated (Um=Umb+Umj) = 3.193760e-01
Total block volume excavated (Vm) = 3.105829e+00
Total change in potential energy (Ub) = 0.000000e+00
Total mass damping work (Wk) = 0.000000e+00
Total viscous boundary work (Wv) = 0.000000e+00
Total friction work (Wj) = 0.000000e+00
Total plastic strain work (Wp) = 0.000000e+00
Total boundary loading work (W) = 3.193760e+01
Total energy stored (Wr=W-Uc-Ub-Wj-Wp) = 8.323573e-01
Total energy released (Wr=Uk+Wk+Wv+Um) = 3.194194e-01
Breakdown of energy stored in joints (Ucj)
------------------------------------------
Total energy stored in tension (Ujt) = 0.000000e+00
Total energy stored in compression (Ujc) = 0.000000e+00
Total energy stored in shear (Ujs) = 0.000000e+00
**********************************************
* Logging ended at Wed Dec 18 22:27:32 2024
**********************************************
References
Salamon, M. D. G. “Energy Considerations in Rock Mechanics: Fundamental Results,” J. S. Afr. Inst. Min. Metall., 84(8), 233-246 (1984).
| Was this helpful? ... | Itasca Software © 2024, Itasca | Updated: Dec 19, 2024 |
