Skip to content

Commit 671104a

Browse files
davidscndavidscn
authored
Support stress data (#125)
Add support for stress data writing
1 parent 42a3d6d commit 671104a

File tree

7 files changed

+246
-6
lines changed

7 files changed

+246
-6
lines changed

‎FSI/FSI.C

+19
Original file line numberDiff line numberDiff line change
@@ -129,6 +129,15 @@ void preciceAdapter::FSI::FluidStructureInteraction::addWriters(std::string data
129129
);
130130
DEBUG(adapterInfo("Added writer: Displacement."));
131131
}
132+
else if(dataName.find("Stress") == 0)
133+
{
134+
interface->addCouplingDataWriter
135+
(
136+
dataName,
137+
new Stress(mesh_, runTime_.timeName(), solverType_) /* TODO: Add any other arguments here */
138+
);
139+
DEBUG(adapterInfo("Added writer: Stress."));
140+
}
132141

133142
// NOTE: If you want to couple another variable, you need
134143
// to add your new coupling data user as a coupling data
@@ -166,6 +175,16 @@ void preciceAdapter::FSI::FluidStructureInteraction::addReaders(std::string data
166175
);
167176
DEBUG(adapterInfo("Added reader: Displacement."));
168177
}
178+
else if(dataName.find("Stress") == 0)
179+
{
180+
interface->addCouplingDataReader
181+
(
182+
dataName,
183+
new Stress(mesh_, runTime_.timeName(), solverType_) /* TODO: Add any other arguments here */
184+
);
185+
DEBUG(adapterInfo("Added reader: Stress."));
186+
}
187+
169188
// NOTE: If you want to couple another variable, you need
170189
// to add your new coupling data user as a coupling data
171190
// writer here (and as a writer above).

‎FSI/FSI.H

+1
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,7 @@
66
#include "FSI/Displacement.H"
77
#include "FSI/DisplacementDelta.H"
88
#include "FSI/Force.H"
9+
#include "FSI/Stress.H"
910

1011
#include "fvCFD.H"
1112

‎FSI/Force.C

+11-2
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,16 @@
22

33
using namespace Foam;
44

5+
preciceAdapter::FSI::Force::Force
6+
(
7+
const Foam::fvMesh& mesh,
8+
const std::string solverType
9+
)
10+
:
11+
mesh_(mesh),solverType_(solverType)
12+
{}
13+
14+
515
preciceAdapter::FSI::Force::Force
616
(
717
const Foam::fvMesh& mesh,
@@ -13,8 +23,7 @@ preciceAdapter::FSI::Force::Force
1323
*/
1424
)
1525
:
16-
mesh_(mesh),
17-
solverType_(solverType)
26+
Force(mesh, solverType)
1827
{
1928
//What about type "basic"?
2029
if (solverType_.compare("incompressible") != 0 && solverType_.compare("compressible") != 0)

‎FSI/Force.H

+10-4
Original file line numberDiff line numberDiff line change
@@ -5,7 +5,6 @@
55

66
#include "fvCFD.H"
77

8-
// TEMPORARY
98
#include "pointFields.H"
109
#include "vectorField.H"
1110
#include "immiscibleIncompressibleTwoPhaseMixture.H"
@@ -31,7 +30,7 @@ private:
3130
//- Force field
3231
Foam::volVectorField * Force_;
3332

34-
33+
protected:
3534
//- Stress tensor (see the OpenFOAM "Forces" function object)
3635
Foam::tmp<Foam::volSymmTensorField> devRhoReff() const;
3736

@@ -56,10 +55,17 @@ public:
5655
*/
5756
);
5857

59-
//- Write the displacement values into the buffer
58+
//- Constructor
59+
Force
60+
(
61+
const Foam::fvMesh& mesh,
62+
const std::string solverType
63+
);
64+
65+
//- Write the forces values into the buffer
6066
void write(double * buffer, bool meshConnectivity, const unsigned int dim);
6167

62-
//- Read the displacement values from the buffer
68+
//- Read the forces values from the buffer
6369
void read(double * buffer, const unsigned int dim);
6470

6571
//- Destructor

‎FSI/Stress.C

+141
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,141 @@
1+
#include "Stress.H"
2+
3+
using namespace Foam;
4+
5+
preciceAdapter::FSI::Stress::Stress
6+
(
7+
const Foam::fvMesh& mesh,
8+
const fileName& timeName,
9+
const std::string solverType
10+
/* TODO: We should add any required field names here.
11+
/ They would need to be vector fields.
12+
/ See FSI/Temperature.C for details.
13+
*/
14+
)
15+
:
16+
Force(mesh,solverType),
17+
mesh_(mesh),
18+
solverType_(solverType)
19+
{
20+
if (solverType_.compare("incompressible") != 0 && solverType_.compare("compressible") != 0)
21+
{
22+
FatalErrorInFunction
23+
<< "Stresses calculation only supports "
24+
<< "compressible or incompressible solver type."
25+
<< exit(FatalError);
26+
}
27+
28+
dataType_ = vector;
29+
30+
Stress_ = new volVectorField
31+
(
32+
IOobject
33+
(
34+
"Stress",
35+
timeName,
36+
mesh,
37+
IOobject::NO_READ,
38+
IOobject::AUTO_WRITE
39+
),
40+
mesh,
41+
dimensionedVector
42+
(
43+
"pdim",
44+
dimensionSet(1,-1,-2,0,0,0,0),
45+
Foam::vector::zero
46+
)
47+
);
48+
}
49+
50+
void preciceAdapter::FSI::Stress::write(double * buffer, bool meshConnectivity, const unsigned int dim)
51+
{
52+
// Compute stress. See the Forces function object.
53+
54+
// Stress tensor boundary field
55+
tmp<volSymmTensorField> tdevRhoReff = this->devRhoReff();
56+
const volSymmTensorField::Boundary& devRhoReffb =
57+
tdevRhoReff().boundaryField();
58+
59+
// Density boundary field
60+
tmp<volScalarField> trho = this->rho();
61+
const volScalarField::Boundary& rhob =
62+
trho().boundaryField();
63+
64+
// Pressure boundary field
65+
tmp<volScalarField> tp = mesh_.lookupObject<volScalarField>("p");
66+
const volScalarField::Boundary& pb =
67+
tp().boundaryField();
68+
69+
int bufferIndex = 0;
70+
// For every boundary patch of the interface
71+
for (uint j = 0; j < patchIDs_.size(); j++)
72+
{
73+
74+
int patchID = patchIDs_.at(j);
75+
76+
// Compute normal vectors on each patch
77+
const vectorField nV = mesh_.boundary()[patchID].nf();
78+
79+
// Pressure constribution
80+
if (solverType_.compare("incompressible") == 0)
81+
{
82+
Stress_->boundaryFieldRef()[patchID] =
83+
nV * pb[patchID] * rhob[patchID];
84+
}
85+
else if (solverType_.compare("compressible") == 0)
86+
{
87+
Stress_->boundaryFieldRef()[patchID] =
88+
nV * pb[patchID];
89+
}
90+
else
91+
{
92+
FatalErrorInFunction
93+
<< "Stress calculation does only support "
94+
<< "compressible or incompressible solver type."
95+
<< exit(FatalError);
96+
}
97+
98+
// Viscous contribution
99+
Stress_->boundaryFieldRef()[patchID] +=
100+
nV & devRhoReffb[patchID];
101+
102+
// Write the stress to the preCICE buffer
103+
// For every cell of the patch
104+
forAll(Stress_->boundaryFieldRef()[patchID], i)
105+
{
106+
// Copy the stress into the buffer
107+
// x-dimension
108+
buffer[bufferIndex++]
109+
=
110+
Stress_->boundaryFieldRef()[patchID][i].x();
111+
112+
// y-dimension
113+
buffer[bufferIndex++]
114+
=
115+
Stress_->boundaryFieldRef()[patchID][i].y();
116+
117+
if(dim == 3)
118+
// z-dimension
119+
buffer[bufferIndex++]
120+
=
121+
Stress_->boundaryFieldRef()[patchID][i].z();
122+
}
123+
}
124+
}
125+
126+
void preciceAdapter::FSI::Stress::read(double * buffer, const unsigned int dim)
127+
{
128+
/* TODO: Implement
129+
* We need two nested for-loops for each patch,
130+
* the outer for the locations and the inner for the dimensions.
131+
* See the preCICE readBlockVectorData() implementation.
132+
*/
133+
FatalErrorInFunction
134+
<< "Reading stresses is not supported."
135+
<< exit(FatalError);
136+
}
137+
138+
preciceAdapter::FSI::Stress::~Stress()
139+
{
140+
delete Stress_;
141+
}

‎FSI/Stress.H

+63
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,63 @@
1+
#ifndef FSI_STRESS_H
2+
#define FSI_STRESS_H
3+
4+
#include "CouplingDataUser.H"
5+
#include "Force.H"
6+
7+
#include "fvCFD.H"
8+
9+
#include "pointFields.H"
10+
#include "vectorField.H"
11+
#include "immiscibleIncompressibleTwoPhaseMixture.H"
12+
#include "turbulentFluidThermoModel.H"
13+
#include "turbulentTransportModel.H"
14+
15+
namespace preciceAdapter
16+
{
17+
namespace FSI
18+
{
19+
20+
//- Class that writes and reads stress [N/m^2]:
21+
// This is essentially a force (in spatial coordinates) scaled by the deformed
22+
// cell face. Thus, a consistent quantity. Calculation concept has been copied
23+
// from the force module, but the scaled version here is commonly used in FEM
24+
// applications.
25+
class Stress : public Force {
26+
27+
private:
28+
29+
//- OpenFOAM fvMesh object (we need to access the objects' registry multiple times)
30+
const Foam::fvMesh& mesh_;
31+
32+
const std::string solverType_;
33+
34+
//- Stress field
35+
Foam::volVectorField * Stress_;
36+
37+
38+
public:
39+
40+
//- Constructor
41+
Stress
42+
(
43+
const Foam::fvMesh& mesh,
44+
const fileName& timeName,
45+
const std::string solverType
46+
// We create an IOobject and we need the time directory
47+
/* TODO: We should add any required field names here.
48+
*/
49+
);
50+
51+
//- Write the stress values into the buffer
52+
void write(double * buffer, bool meshConnectivity, const unsigned int dim);
53+
54+
//- Read the stress values from the buffer
55+
void read(double * buffer, const unsigned int dim);
56+
57+
//- Destructor
58+
~Stress();
59+
};
60+
}
61+
}
62+
63+
#endif

‎Make/files

+1
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,7 @@ CHT/CHT.C
1414

1515
FSI/FSI.C
1616
FSI/Force.C
17+
FSI/Stress.C
1718
FSI/Displacement.C
1819
FSI/DisplacementDelta.C
1920

0 commit comments

Comments
 (0)