Line data Source code
1 : /*
2 : * Copyright (c) 2013-2015: G-CSC, Goethe University Frankfurt
3 : * Author: Raphael Prohl
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 : #ifndef __H__UG__LIB_ALGEBRA__OPERATOR__PRECONDITIONER__PROJECTED_GAUSS_SEIDEL__SCALAR_OBSTACLE_IMPL__
34 : #define __H__UG__LIB_ALGEBRA__OPERATOR__PRECONDITIONER__PROJECTED_GAUSS_SEIDEL__SCALAR_OBSTACLE_IMPL__
35 :
36 : #include "scalar_obstacle.h"
37 :
38 : namespace ug{
39 :
40 : //////////////////////////////
41 : // SCALAR LOWER OBSTACLE
42 : //////////////////////////////
43 :
44 : template <typename TDomain, typename TAlgebra>
45 : void
46 0 : ScalarLowerObstacle<TDomain, TAlgebra>::
47 : adjust_sol_and_cor(value_type& sol_i, value_type& c_i, bool& dofIsActive,
48 : const DoFIndex& dof)
49 : {
50 0 : const size_t comp = dof[1];
51 :
52 : // tmpSol := u_{s-1/2} = u_{s-1} + c
53 0 : const number tmpSol = BlockRef(sol_i, comp) + BlockRef(c_i, comp);
54 :
55 : // get lower obstacle value corresponding to the dof
56 0 : const number obsVal = m_mObstacleValues[dof];
57 :
58 : // check, if dof is active (tmpSol <= obsVal)
59 0 : if (!(tmpSol > obsVal))
60 : {
61 : // is active DoF
62 0 : m_vActiveDofs.push_back(dof);
63 :
64 : // adjust correction & set solution to obstacle-value
65 0 : BlockRef(c_i, comp) = obsVal - BlockRef(sol_i, comp);
66 0 : BlockRef(sol_i, comp) = obsVal;
67 0 : dofIsActive = true;
68 : }
69 0 : }
70 :
71 : template <typename TDomain, typename TAlgebra>
72 : void
73 0 : ScalarLowerObstacle<TDomain, TAlgebra>::
74 : adjust_defect_to_constraint(vector_type& d)
75 : {
76 0 : for (std::vector<MultiIndex<2> >::iterator itActiveInd = m_vActiveDofs.begin();
77 0 : itActiveInd < m_vActiveDofs.end(); ++itActiveInd)
78 : {
79 : // check, if Ax <= b. For that case the new defect is set to zero,
80 : // since all equations/constraints are fulfilled
81 0 : number defect = BlockRef(d[(*itActiveInd)[0]], (*itActiveInd)[1]);
82 0 : if (defect < 0.0)
83 0 : BlockRef(d[(*itActiveInd)[0]], (*itActiveInd)[1]) = 0.0;
84 : }
85 0 : }
86 :
87 : template <typename TDomain, typename TAlgebra>
88 : void
89 0 : ScalarLowerObstacle<TDomain, TAlgebra>::
90 : restrict_obs_values()
91 0 : {}
92 :
93 : //////////////////////////////
94 : // SCALAR UPPER OBSTACLE
95 : //////////////////////////////
96 :
97 : template <typename TDomain, typename TAlgebra>
98 : void
99 0 : ScalarUpperObstacle<TDomain, TAlgebra>::
100 : adjust_sol_and_cor(value_type& sol_i, value_type& c_i, bool& dofIsActive,
101 : const DoFIndex& dof)
102 : {
103 0 : const size_t comp = dof[1];
104 :
105 : // tmpSol := u_{s-1/2} = u_{s-1} + c
106 0 : const number tmpSol = BlockRef(sol_i, comp) + BlockRef(c_i, comp);
107 :
108 : // get upper obstacle value corresponding to the dof
109 0 : const number obsVal = m_mObstacleValues[dof];
110 :
111 : // check, if dof is active (tmpSol >= obsVal)
112 0 : if (!(tmpSol < obsVal))
113 : {
114 : // is active DoF
115 0 : m_vActiveDofs.push_back(dof);
116 :
117 : // adjust correction & set solution to obstacle-value
118 0 : BlockRef(c_i, comp) = obsVal - BlockRef(sol_i, comp);
119 0 : BlockRef(sol_i, comp) = obsVal;
120 0 : dofIsActive = true;
121 : }
122 : //UG_LOG("dof " <<dof<< " is active: " <<dofIsActive<<"\n");
123 0 : }
124 :
125 : template <typename TDomain, typename TAlgebra>
126 : void
127 0 : ScalarUpperObstacle<TDomain, TAlgebra>::
128 : adjust_defect_to_constraint(vector_type& d)
129 : {
130 0 : for (std::vector<MultiIndex<2> >::iterator itActiveInd = m_vActiveDofs.begin();
131 0 : itActiveInd < m_vActiveDofs.end(); ++itActiveInd)
132 : {
133 : // check, if Ax > b. For that case the new defect is set to zero,
134 : // since all equations/constraints are fulfilled
135 : //UG_LOG("adjust_defect: " << (*itActiveInd)[0] <<","<< (*itActiveInd)[1] << "\n");
136 0 : number defect = BlockRef(d[(*itActiveInd)[0]], (*itActiveInd)[1]);
137 0 : if (defect > 0.0)
138 : {
139 : //UG_LOG("defect > 0 \n");
140 0 : BlockRef(d[(*itActiveInd)[0]], (*itActiveInd)[1]) = 0.0;
141 : }
142 : }
143 0 : }
144 :
145 : template <typename TDomain, typename TAlgebra>
146 : void
147 0 : ScalarUpperObstacle<TDomain, TAlgebra>::
148 : restrict_obs_values()
149 0 : {}
150 :
151 : } // end namespace ug
152 :
153 : #endif /* __H__UG__LIB_ALGEBRA__OPERATOR__PRECONDITIONER__PROJECTED_GAUSS_SEIDEL__SCALAR_OBSTACLE_IMPL__ */
|
