Skip to content

Commit 23b1baf

Browse files
datafl4shdatafl4sh
committed
New code for obstacle problem
1 parent baa7316 commit 23b1baf

16 files changed

Lines changed: 1337 additions & 462 deletions

File tree

‎apps/CMakeLists.txt‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -15,7 +15,7 @@ add_subdirectory(linear_elasticity)
1515
#add_subdirectory(fvca5_meshinfo)
1616
add_subdirectory(integration)
1717
#add_subdirectory(mesh_new)
18-
#add_subdirectory(obstacle)
18+
add_subdirectory(obstacle)
1919
add_subdirectory(contact)
2020
#add_subdirectory(test_bases)
2121
add_subdirectory(new_code_test)

‎apps/bingham/bingham_vector_solver.hpp‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -91,7 +91,7 @@ class ADDM
9191
auto value = 1./(2. * (vp.mu + vp.alpha));
9292
auto G = disk::make_hlow_stokes(msh, cl, di, true);
9393

94-
auto cell_ofs = disk::priv::offset(msh, cl);
94+
auto cell_ofs = disk::offset(msh, cl);
9595
vector_type Gu = G.first * u_TF;
9696
vector_type stress = multiplier.block(cell_ofs * sbs, 0, sbs, 1);
9797

@@ -127,7 +127,7 @@ class ADDM
127127
auto cb = disk::make_vector_monomial_basis(msh, cl, di.cell_degree());
128128
auto sb = disk::make_sym_matrix_monomial_basis(msh, cl, di.face_degree());
129129

130-
auto cell_ofs = disk::priv::offset(msh, cl);
130+
auto cell_ofs = disk::offset(msh, cl);
131131
auto num_faces = howmany_faces(msh, cl);
132132

133133
vector_type rhs = vector_type::Zero(cbs + fbs * num_faces);

‎apps/bingham/viscoplasticity_utils.hpp‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -425,7 +425,7 @@ plot_over_line(const Mesh & msh,
425425
if(wn_PnPoly( msh, cl, p))
426426
{
427427
auto cbs = disk::vector_basis_size(cell_degree, dim, dim);
428-
auto cell_ofs = disk::priv::offset(msh, cl);
428+
auto cell_ofs = disk::offset(msh, cl);
429429

430430
vector_type s = vec.block(cell_ofs * cbs, 0, cbs, 1);
431431
auto cb = disk::make_vector_monomial_basis(msh, cl, cell_degree);
@@ -464,7 +464,7 @@ find_values_at_points(const Mesh & msh,
464464
if(wn_PnPoly( msh, cl, p))
465465
{
466466
auto cbs = vector_basis_size(cell_degree, dim, dim);
467-
auto cell_ofs = disk::priv::offset(msh, cl);
467+
auto cell_ofs = disk::offset(msh, cl);
468468

469469
vector_type s = vec.block(cell_ofs * cbs, 0, cbs, 1);
470470
auto cb = make_vector_monomial_basis(msh, cl, cell_degree);
@@ -510,7 +510,7 @@ compute_discontinuous_velocity(const Mesh& msh,
510510

511511
for (auto& cl : msh)
512512
{
513-
auto cell_ofs = disk::priv::offset(msh, cl);
513+
auto cell_ofs = disk::offset(msh, cl);
514514
Matrix<T, Dynamic, 1> x = sol.block(cell_ofs * cbs, 0, cbs, 1);
515515

516516
const auto cell_nodes = post_mesh.nodes_cell(cell_i);
@@ -573,7 +573,7 @@ quiver( const Mesh& msh, const dynamic_vector<T>& sol, const Assembler& assemble
573573
Mesh::dimension, Mesh::dimension);
574574
for (auto& cl: msh)
575575
{
576-
auto cell_ofs = disk::priv::offset(msh, cl);
576+
auto cell_ofs = disk::offset(msh, cl);
577577
Matrix<T, Dynamic, 1> s = assembler.take_velocity(msh, cl, sol);
578578

579579
auto cb = disk::make_vector_monomial_basis(msh, cl, di.cell_degree());

‎apps/contact/common.hpp‎

Lines changed: 22 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -1881,7 +1881,7 @@ class diffusion_condensed_assembler_nitsche_cells
18811881
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
18821882
{
18831883
auto fc = fcs[face_i];
1884-
auto face_offset = priv::offset(msh, fc);
1884+
auto face_offset = offset(msh, fc);
18851885
auto face_LHS_offset = compress_table.at(face_offset)*fbs;
18861886

18871887
auto face_id = msh.lookup(fc);
@@ -1942,7 +1942,7 @@ class diffusion_condensed_assembler_nitsche_cells
19421942
}
19431943
else
19441944
{
1945-
auto face_offset = priv::offset(msh, fc);
1945+
auto face_offset = offset(msh, fc);
19461946
auto face_SOL_offset = compress_table.at(face_offset)*fbs;
19471947
ret.block(face_i*fbs, 0, fbs, 1) = solution.block(face_SOL_offset, 0, fbs, 1);
19481948
}
@@ -2076,7 +2076,7 @@ class diffusion_condensed_assembler_nitsche_faces
20762076
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
20772077
{
20782078
auto fc = fcs[face_i];
2079-
auto face_offset = priv::offset(msh, fc);
2079+
auto face_offset = offset(msh, fc);
20802080
auto face_LHS_offset = compress_table.at(face_offset)*fbs;
20812081

20822082
auto face_id = msh.lookup(fc);
@@ -2125,7 +2125,7 @@ class diffusion_condensed_assembler_nitsche_faces
21252125
if (!eid.first) throw std::invalid_argument("This is a bug: face not found");
21262126
const auto face_id=eid.second;
21272127

2128-
auto face_offset = priv::offset(msh, fc);
2128+
auto face_offset = offset(msh, fc);
21292129
auto face_SOL_offset = compress_table.at(face_offset)*fbs;
21302130
ret.block(face_i*fbs, 0, fbs, 1) = solution.block(face_SOL_offset, 0, fbs, 1);
21312131
}
@@ -2263,7 +2263,7 @@ class diffusion_full_assembler
22632263
std::vector<assembly_index> asm_map;
22642264
asm_map.reserve(cbs + fcs.size()*fbs);
22652265

2266-
auto cell_offset = priv::offset(msh, cl);
2266+
auto cell_offset = offset(msh, cl);
22672267
auto cell_LHS_offset = cell_offset * cbs;
22682268

22692269
for (size_t i = 0; i < cbs; i++)
@@ -2272,7 +2272,7 @@ class diffusion_full_assembler
22722272
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
22732273
{
22742274
auto fc = fcs[face_i];
2275-
auto face_offset = priv::offset(msh, fc);
2275+
auto face_offset = offset(msh, fc);
22762276
auto face_LHS_offset = cbs * msh.cells_size() + compress_table.at(face_offset)*fbs;
22772277

22782278
auto fc_id = msh.lookup(fc);
@@ -2328,7 +2328,7 @@ class diffusion_full_assembler
23282328
{
23292329
auto num_faces = howmany_faces(msh, cl);
23302330
auto dim = Mesh::dimension;
2331-
auto cell_ofs = disk::priv::offset(msh, cl);
2331+
auto cell_ofs = disk::offset(msh, cl);
23322332

23332333
Matrix<T, Dynamic, 1> svel(cbs + num_faces * fbs );
23342334
svel.block(0, 0, cbs, 1) = sol.block(cell_ofs * cbs, 0, cbs, 1);
@@ -2350,7 +2350,7 @@ class diffusion_full_assembler
23502350
}
23512351
else
23522352
{
2353-
auto face_ofs = priv::offset(msh, fc);
2353+
auto face_ofs = offset(msh, fc);
23542354
auto global_ofs = cbs * msh.cells_size() + compress_table.at(face_ofs)*fbs;
23552355
svel.block(cbs + i*fbs, 0, fbs, 1) = sol.block(global_ofs, 0, fbs, 1);
23562356
}
@@ -2508,7 +2508,7 @@ class diffusion_mix_full_assembler
25082508
std::vector<assembly_index> asm_map;
25092509
asm_map.reserve(cbs + fcs.size()*fbs);
25102510

2511-
auto cell_offset = priv::offset(msh, cl);
2511+
auto cell_offset = offset(msh, cl);
25122512
auto cell_LHS_offset = cell_offset * cbs;
25132513

25142514
element_type ecell = CELL;
@@ -2521,7 +2521,7 @@ class diffusion_mix_full_assembler
25212521
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
25222522
{
25232523
auto fc = fcs[face_i];
2524-
auto face_offset = priv::offset(msh, fc);
2524+
auto face_offset = offset(msh, fc);
25252525
auto face_LHS_offset = cbs * msh.cells_size() + compress_table.at(face_offset)*fbs;
25262526

25272527
auto fc_id = msh.lookup(fc);
@@ -2648,7 +2648,7 @@ class diffusion_mix_full_assembler
26482648
{
26492649
auto num_faces = howmany_faces(msh, cl);
26502650
auto dim = Mesh::dimension;
2651-
auto cell_ofs = disk::priv::offset(msh, cl);
2651+
auto cell_ofs = disk::offset(msh, cl);
26522652

26532653
Matrix<T, Dynamic, 1> svel(cbs + num_faces * fbs );
26542654
svel.block(0, 0, cbs, 1) = sol.block(cell_ofs * cbs, 0, cbs, 1);
@@ -2672,7 +2672,7 @@ class diffusion_mix_full_assembler
26722672
svel.block(cbs + i * fbs, 0, fbs, 1) = Matrix<T, Dynamic, 1>::Zero(fbs);
26732673
else
26742674
{
2675-
auto face_ofs = priv::offset(msh, fc);
2675+
auto face_ofs = offset(msh, fc);
26762676
auto global_ofs = cbs * msh.cells_size() + compress_table.at(face_ofs)*fbs;
26772677
svel.block(cbs + i*fbs, 0, fbs, 1) = sol.block(global_ofs, 0, fbs, 1);
26782678
}
@@ -2833,7 +2833,7 @@ class contact_full_assembler
28332833
std::vector<assembly_index> asm_map;
28342834
asm_map.reserve(cbs + fcs.size()*fbs);
28352835

2836-
auto cell_offset = priv::offset(msh, cl);
2836+
auto cell_offset = offset(msh, cl);
28372837
auto cell_LHS_offset = cell_offset * cbs;
28382838

28392839
element_type ecell = CELL;
@@ -2845,7 +2845,7 @@ class contact_full_assembler
28452845
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
28462846
{
28472847
auto fc = fcs[face_i];
2848-
auto face_offset = priv::offset(msh, fc);
2848+
auto face_offset = offset(msh, fc);
28492849
auto face_LHS_offset = cbs * msh.cells_size() + compress_table.at(face_offset)*fbs;
28502850

28512851
auto fc_id = msh.lookup(fc);
@@ -2969,7 +2969,7 @@ class contact_full_assembler
29692969
{
29702970
auto num_faces = howmany_faces(msh, cl);
29712971
auto dim = Mesh::dimension;
2972-
auto cell_ofs = disk::priv::offset(msh, cl);
2972+
auto cell_ofs = disk::offset(msh, cl);
29732973
auto num_total_dofs = cbs + num_faces * fbs;
29742974
Matrix<T, Dynamic, 1> svel = Matrix<T, Dynamic, 1>::Zero(num_total_dofs);
29752975
svel.block(0, 0, cbs, 1) = sol.block(cell_ofs * cbs, 0, cbs, 1);
@@ -2996,7 +2996,7 @@ class contact_full_assembler
29962996
svel.block(cbs + i * fbs, 0, fbs, 1) = vector_type::Zero(fbs);
29972997
else
29982998
{
2999-
auto face_ofs = priv::offset(msh, fc);
2999+
auto face_ofs = offset(msh, fc);
30003000
auto global_ofs = cbs * msh.cells_size() + compress_table.at(face_ofs)*fbs;
30013001
svel.block(cbs + i*fbs, 0, fbs, 1) = sol.block(global_ofs, 0, fbs, 1);
30023002
}
@@ -3150,7 +3150,7 @@ class contact_full_assembler_new
31503150
std::vector<assembly_index> asm_map;
31513151
asm_map.reserve(cbs + fcs.size()*fbs);
31523152

3153-
auto cell_offset = priv::offset(msh, cl);
3153+
auto cell_offset = offset(msh, cl);
31543154
auto cell_LHS_offset = cell_offset * cbs;
31553155

31563156
element_type ecell = CELL;
@@ -3160,7 +3160,7 @@ class contact_full_assembler_new
31603160
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
31613161
{
31623162
auto fc = fcs[face_i];
3163-
auto face_offset = priv::offset(msh, fc);
3163+
auto face_offset = offset(msh, fc);
31643164
auto face_LHS_offset = cbs * msh.cells_size() + compress_table.at(face_offset)*fbs;
31653165

31663166
auto fc_id = msh.lookup(fc);
@@ -3293,7 +3293,7 @@ class contact_full_assembler_new
32933293
{
32943294
auto num_faces = howmany_faces(msh, cl);
32953295
auto dim = Mesh::dimension;
3296-
auto cell_ofs = disk::priv::offset(msh, cl);
3296+
auto cell_ofs = disk::offset(msh, cl);
32973297
auto num_total_dofs = cbs + num_faces * fbs;
32983298
Matrix<T, Dynamic, 1> svel = Matrix<T, Dynamic, 1>::Zero(num_total_dofs);
32993299
svel.block(0, 0, cbs, 1) = sol.block(cell_ofs * cbs, 0, cbs, 1);
@@ -3311,7 +3311,7 @@ class contact_full_assembler_new
33113311
svel.block(cbs + i * fbs, 0, fbs, 1) = Matrix<T, Dynamic, 1>::Zero(fbs);
33123312
else
33133313
{
3314-
auto face_ofs = priv::offset(msh, fc);
3314+
auto face_ofs = offset(msh, fc);
33153315
auto global_ofs = cbs * msh.cells_size() + compress_table.at(face_ofs)*fbs;
33163316
svel.block(cbs + i*fbs, 0, fbs, 1) = sol.block(global_ofs, 0, fbs, 1);
33173317
}
@@ -3436,7 +3436,7 @@ class contact_face_assembler_new
34363436
for (size_t face_i = 0; face_i < fcs.size(); face_i++)
34373437
{
34383438
auto fc = fcs[face_i];
3439-
auto face_offset = priv::offset(msh, fc);
3439+
auto face_offset = offset(msh, fc);
34403440
auto face_LHS_offset = compress_table.at(face_offset)*fbs;
34413441

34423442
auto face_id = msh.lookup(fc);
@@ -3616,7 +3616,7 @@ class contact_face_assembler_new
36163616
ret.block(face_i * fbs, 0, fbs, 1) = Matrix<T, Dynamic, 1>::Zero(fbs);
36173617
else
36183618
{
3619-
auto face_offset = priv::offset(msh, fc);
3619+
auto face_offset = offset(msh, fc);
36203620
auto face_SOL_offset = compress_table.at(face_offset)*fbs;
36213621
ret.block(face_i * fbs, 0, fbs, 1) = solution.block(face_SOL_offset, 0, fbs, 1);
36223622
}

‎apps/diffusion/diffusion_alg.cpp‎

Lines changed: 7 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -129,7 +129,7 @@ class augmented_lagrangian_diffusion
129129
{
130130
auto bar = barycenter(msh, cl);
131131
Matrix<scalar_type, Dynamic, 1> p = project_function(msh, cl, di, sol_fun);
132-
auto cell_ofs = disk::priv::offset(msh, cl);
132+
auto cell_ofs = disk::offset(msh, cl);
133133
Matrix<scalar_type, Dynamic, 1> s = sol.block(cell_ofs * cbs, 0, cbs, 1);
134134
Matrix<scalar_type, Dynamic, 1> diff = s - p.head(cbs);
135135
auto cb = disk::make_scalar_monomial_basis(msh, cl, di.cell_degree());
@@ -161,7 +161,7 @@ class augmented_lagrangian_diffusion
161161
const cell_type& cl)
162162
{
163163
auto gamma = compute_auxiliar(msh, cl);
164-
auto cell_ofs = disk::priv::offset(msh, cl);
164+
auto cell_ofs = disk::offset(msh, cl);
165165
auxiliar.block(cell_ofs * sbs, 0, sbs, 1) = gamma;
166166
return;
167167
}
@@ -172,7 +172,7 @@ class augmented_lagrangian_diffusion
172172
{
173173
auto sb = disk::make_vector_monomial_basis(msh, cl, di.cell_degree());
174174

175-
auto cell_ofs = disk::priv::offset(msh, cl);
175+
auto cell_ofs = disk::offset(msh, cl);
176176
auto num_total_dofs = cbs + howmany_faces(msh, cl) * fbs;
177177

178178
Eigen::Matrix<scalar_type, Eigen::Dynamic, 1> u_TF =
@@ -196,7 +196,7 @@ class augmented_lagrangian_diffusion
196196

197197
for(auto cl: msh)
198198
{
199-
auto cell_ofs = disk::priv::offset(msh, cl);
199+
auto cell_ofs = disk::offset(msh, cl);
200200
auto num_total_dofs = cbs + howmany_faces(msh, cl) * fbs;
201201

202202
Eigen::Matrix<scalar_type, Eigen::Dynamic, 1> u_TF =
@@ -228,7 +228,7 @@ class augmented_lagrangian_diffusion
228228
auto G = make_vector_hho_gradrec(msh, cl, di);
229229
auto cb = disk::make_scalar_monomial_basis(msh, cl, di.cell_degree());
230230
auto sb = disk::make_vector_monomial_basis(msh, cl, di.cell_degree());
231-
auto cell_ofs = disk::priv::offset(msh, cl);
231+
auto cell_ofs = disk::offset(msh, cl);
232232
auto num_faces = howmany_faces(msh, cl);
233233

234234
auto stress = multiplier.block( sbs * cell_ofs, 0, sbs, 1);
@@ -300,7 +300,7 @@ class augmented_lagrangian_diffusion
300300
auto gr = make_scalar_hho_laplacian(msh, cl, di);
301301
auto stab = make_scalar_hho_stabilization(msh, cl, gr.first, di);
302302

303-
auto cell_ofs = disk::priv::offset(msh, cl);
303+
auto cell_ofs = disk::offset(msh, cl);
304304
auto num_total_dofs = cbs + fbs * howmany_faces(msh, cl);
305305

306306
Matrix<scalar_type, Dynamic, 1> cell_rhs = rhs_all.block( cbs*cell_ofs, 0, cbs, 1);
@@ -332,7 +332,7 @@ class augmented_lagrangian_diffusion
332332

333333
for (auto& cl: msh)
334334
{
335-
auto cell_ofs = disk::priv::offset(msh, cl);
335+
auto cell_ofs = disk::offset(msh, cl);
336336
vector_type stress = multiplier.block(cell_ofs * sbs, 0, sbs, 1);
337337
auto gamma = compute_auxiliar( msh, cl);
338338

‎apps/obstacle/CMakeLists.txt‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,2 +1,6 @@
1+
option(BUILD_OBSTACLE "Build Obstacle Problem application" ON)
2+
if (BUILD_OBSTACLE)
13
add_executable(obstacle obstacle.cpp)
24
target_link_libraries(obstacle ${LINK_LIBS})
5+
endif (BUILD_OBSTACLE)
6+

‎apps/obstacle/gnuplot_scripts/nitsche_2d.plot‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
set term postscript enhanced color font 'Helvetica,22'
2+
set output 'nitsche_2d.eps'
3+
4+
set logscale x
5+
set logscale y
6+
7+
set grid xtics mxtics ytics
8+
9+
set xrange [1:0.01] reverse
10+
11+
set style line 1 lt 1 lc rgb "#006600" lw 1
12+
set style line 2 lt 3 dt 2 lc rgb "#990000" lw 1
13+
14+
plot 'triangles_k0.txt' using 1:4 w lp ls 1 pt 8 ps 1.5 title "Triangles, k=0", \
15+
'quads_k0.txt' using 1:4 w lp ls 1 pt 4 ps 1.5 title "Quads, k=0", \
16+
'hexagons_k0.txt' using 1:4 w lp ls 1 pt 14 ps 1.5 title "Hexagons, k=0", \
17+
'triangles_k1.txt' using 1:4 w lp ls 2 pt 8 ps 1.5 title "Triangles, k=1", \
18+
'quads_k1.txt' using 1:4 w lp ls 2 pt 4 ps 1.5 title "Quads, k=1", \
19+
'hexagons_k1.txt' using 1:4 w lp ls 2 pt 14 ps 1.5 title "Hexagons, k=1"

‎apps/obstacle/gnuplot_scripts/nitsche_3d.plot‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
set term postscript enhanced color font 'Helvetica,22'
2+
set output 'nitsche_3d.eps'
3+
4+
set logscale x
5+
set logscale y
6+
7+
set grid xtics mxtics ytics
8+
9+
set xrange [1:0.01] reverse
10+
11+
set style line 1 lt 1 lc rgb "#006600" lw 1
12+
set style line 2 lt 3 dt 2 lc rgb "#990000" lw 1
13+
14+
plot 'tetras_k0.txt' using 1:4 w lp ls 1 pt 8 ps 1.5 title "Tetrahedrons, k=0", \
15+
'hexas_k0.txt' using 1:4 w lp ls 1 pt 4 ps 1.5 title "Cubes, k=0", \
16+
'prisms_k0.txt' using 1:4 w lp ls 1 pt 14 ps 1.5 title "Prisms, k=0", \
17+
'tetras_k1.txt' using 1:4 w lp ls 2 pt 8 ps 1.5 title "Tetrahedrons, k=1", \
18+
'hexas_k1.txt' using 1:4 w lp ls 2 pt 4 ps 1.5 title "Cubes, k=1", \
19+
'prisms_k1.txt' using 1:4 w lp ls 2 pt 14 ps 1.5 title "Prisms, k=1"

‎apps/obstacle/gnuplot_scripts/strong_2d.plot‎

Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
set term postscript enhanced color font 'Helvetica,22'
2+
set output 'strong_2d.eps'
3+
4+
set logscale x
5+
set logscale y
6+
7+
set grid xtics mxtics ytics
8+
9+
set xrange [1:0.01] reverse
10+
11+
set style line 1 lt 1 lc rgb "#006600" lw 1
12+
set style line 2 lt 3 dt 2 lc rgb "#990000" lw 1
13+
14+
plot 'triangles_k0.txt' using 1:2 w lp ls 1 pt 8 ps 1.5 title "Triangles, k=0", \
15+
'quads_k0.txt' using 1:2 w lp ls 1 pt 4 ps 1.5 title "Quads, k=0", \
16+
'hexagons_k0.txt' using 1:2 w lp ls 1 pt 14 ps 1.5 title "Hexagons, k=0", \
17+
'triangles_k1.txt' using 1:2 w lp ls 2 pt 8 ps 1.5 title "Triangles, k=1", \
18+
'quads_k1.txt' using 1:2 w lp ls 2 pt 4 ps 1.5 title "Quads, k=1", \
19+
'hexagons_k1.txt' using 1:2 w lp ls 2 pt 14 ps 1.5 title "Hexagons, k=1"

0 commit comments

Comments
 (0)