Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F81847913
parabolic_driver.hpp
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sun, Sep 8, 06:11
Size
6 KB
Mime Type
text/x-c++
Expires
Tue, Sep 10, 06:11 (1 d, 21 h)
Engine
blob
Format
Raw Data
Handle
20626242
Attached To
rSPECMICP SpecMiCP / ReactMiCP
parabolic_driver.hpp
View Options
/* =============================================================================
Copyright (c) 2014 - 2016
F. Georget <fabieng@princeton.edu> Princeton University
All rights reserved.
Redistribution and use in source and binary forms, with or without modification,
are permitted provided that the following conditions are met:
1. Redistributions of source code must retain the above copyright notice,
this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
3. Neither the name of the copyright holder nor the names of its
contributors may be used to endorse or promote products derived from this
software without specific prior written permission.
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR
ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
============================================================================= */
#ifndef SPECMICP_DFPMSOLVER_PARABOLICDRIVER_HPP
#define SPECMICP_DFPMSOLVER_PARABOLICDRIVER_HPP
#include "specmicp_common/types.hpp"
#include "specmicp_common/sparse_solvers/sparse_solver.hpp"
#include "driver.hpp"
#include "parabolic_structs.hpp"
namespace
specmicp
{
namespace
dfpmsolver
{
//! \brief The parabolic driver
//!
//! Parabolic driver for finite element programs
template
<
class
Program
>
class
ParabolicDriver
:
public
Driver
<
Program
,
ParabolicDriverOptions
,
ParabolicDriverPerformance
>
{
public
:
using
base
=
Driver
<
Program
,
ParabolicDriverOptions
,
ParabolicDriverPerformance
>
;
using
base
::
program
;
using
base
::
get_neq
;
using
base
::
get_options
;
using
base
::
get_perfs
;
using
base
::
scaling
;
using
base
::
initialize_scaling
;
ParabolicDriver
(
Program
&
the_program
)
:
Driver
<
Program
,
ParabolicDriverOptions
,
ParabolicDriverPerformance
>
(
the_program
),
m_velocity
(
Vector
::
Zero
(
the_program
.
get_tot_ndf
())),
m_solver
(
nullptr
)
{}
//! \brief Solve a timestep of length dt
ParabolicDriverReturnCode
solve_timestep
(
scalar_t
dt
,
Vector
&
displacement
);
//! \brief Restart the current timestep
ParabolicDriverReturnCode
restart_timestep
(
Vector
&
displacement
);
//! \brief Check if the solution has been found
ParabolicDriverReturnCode
check_convergence
();
//! \brief Compute the residuals
void
compute_residuals
(
const
Vector
&
displacement
,
Vector
&
residual
)
{
compute_residuals
(
displacement
,
m_velocity
,
residual
);
}
//! \brief Compute the residuals
void
compute_residuals
(
const
Vector
&
displacement
,
const
Vector
&
velocity
,
Vector
&
residual
)
{
program
().
compute_residuals
(
displacement
,
velocity
,
residual
);
get_perfs
().
nb_call_residuals
+=
1
;
}
//! \brief Compute the residuals for use to check convergence
//!
//! For most program, it is the same as calling compute_residuals
void
compute_residuals_0
(
const
Vector
&
displacement
,
const
Vector
&
velocity
,
Vector
&
residual
)
{
program
().
compute_residuals_0
(
displacement
,
velocity
,
residual
);
get_perfs
().
nb_call_residuals
+=
1
;
}
//! \brief Compute the jacobian
void
compute_jacobian
(
Vector
&
displacement
,
Vector
&
velocity
,
Eigen
::
SparseMatrix
<
scalar_t
>&
jacobian
);
//! \brief Return the norm of the current residuals
double
norm_residuals
()
{
return
m_residuals
.
norm
();}
//! \brief Read/Write reference to the velocity vector
const
Vector
&
get_velocity
()
const
{
return
m_velocity
;}
Vector
&
velocity
()
{
return
m_velocity
;}
void
set_velocity
(
Vector
&
velocity_vector
);
//! \brief Reset the Sparse solver
//! This function should be called when a new solver need to be used
void
reset_solver
()
{
return
m_solver
.
reset
(
nullptr
);}
//! \brief Call this function if the pattern of the jacobian has changed
void
jacobian_pattern_has_changed
()
{
reset_solver
();}
//! \brief Initialize the computation
void
initialize_timestep
(
scalar_t
dt
,
Eigen
::
VectorXd
&
displacement
);
//! \brief Strang Linesearch
//!
//! ref :
//! - Matthies et al. (1979)
//! - JHP course notes
ParabolicLinesearchReturnCode
strang_linesearch
(
Vector
&
update
,
Vector
&
displacements
,
scalar_t
&
lambda
);
scalar_t
residuals_0
()
{
return
m_norm_0
;}
private
:
void
reset_velocity
();
//! \brief Backtracking Linesearch
//!
//! ref :
//! - Algo A6.3.1 : Dennis and Schnabel (1983)
//! - Nocedal & Wrigth (2006)
ParabolicLinesearchReturnCode
backtracking_linesearch
(
Vector
&
update
,
Vector
&
displacements
,
scalar_t
&
lambda_out
);
//! Update the variables in displacement
void
update_variable
(
const
Vector
&
update
,
scalar_t
lambda
,
Vector
&
displacement
);
//! \brief Set the predictor
void
set_predictor
(
Vector
&
displacement
);
//! \brief perform the linesearch
ParabolicDriverReturnCode
linesearch
(
Vector
&
update
,
Vector
&
displacements
);
//! Compute the residuals for the linesearch
double
compute_residuals_linesearch
(
Vector
&
update
,
scalar_t
lambda
,
Vector
&
displacement
);
double
compute_residuals_strang_linesearch
(
Vector
&
update
,
scalar_t
lambda
,
Vector
&
displacement
);
scalar_t
update_norm
(
const
Vector
&
update
);
double
m_norm_0
{
0.0
};
double
m_current_dt
{
-
1.0
};
Eigen
::
VectorXd
m_gradient
;
Eigen
::
VectorXd
m_velocity
;
Eigen
::
VectorXd
m_residuals
;
Eigen
::
VectorXd
m_predictor
;
Eigen
::
SparseMatrix
<
double
>
m_jacobian
;
sparse_solvers
::
SparseSolverPtr
<
Eigen
::
SparseMatrix
<
double
>
,
Vector
,
Vector
>
m_solver
;
bool
m_velocity_is_initialized
{
false
};
};
}
// end namespace dfpmsolver
}
// end namespace specmicp
// implementation
// ==============
#include "parabolic_driver.inl"
#endif
// SPECMICP_DFPMSOLVER_PARABOLICDRIVER_HPP
Event Timeline
Log In to Comment