lensproblem.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
28 #ifndef EWOMS_LENS_PROBLEM_HH
29 #define EWOMS_LENS_PROBLEM_HH
30 
34 
35 #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
36 #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
37 #include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
38 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
39 #include <opm/material/fluidsystems/TwoPhaseImmiscibleFluidSystem.hpp>
40 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
41 #include <opm/material/components/SimpleH2O.hpp>
42 #include <opm/material/components/Dnapl.hpp>
43 #include <opm/common/Unused.hpp>
44 
45 #include <dune/common/version.hh>
46 #include <dune/common/fvector.hh>
47 #include <dune/common/fmatrix.hh>
48 
49 #include <sstream>
50 #include <string>
51 #include <iostream>
52 
53 namespace Ewoms {
54 template <class TypeTag>
56 
57 namespace Properties {
59 
60 // declare the properties specific for the lens problem
61 NEW_PROP_TAG(LensLowerLeftX);
62 NEW_PROP_TAG(LensLowerLeftY);
63 NEW_PROP_TAG(LensLowerLeftZ);
64 NEW_PROP_TAG(LensUpperRightX);
65 NEW_PROP_TAG(LensUpperRightY);
66 NEW_PROP_TAG(LensUpperRightZ);
67 
68 // Set the problem property
69 SET_TYPE_PROP(LensBaseProblem, Problem, Ewoms::LensProblem<TypeTag>);
70 
71 // Use Dune-grid's YaspGrid
72 SET_TYPE_PROP(LensBaseProblem, Grid, Dune::YaspGrid<2>);
73 
74 // Set the wetting phase
75 SET_PROP(LensBaseProblem, WettingPhase)
76 {
77 private:
78  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
79 
80 public:
81  typedef Opm::LiquidPhase<Scalar, Opm::SimpleH2O<Scalar> > type;
82 };
83 
84 // Set the non-wetting phase
85 SET_PROP(LensBaseProblem, NonwettingPhase)
86 {
87 private:
88  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
89 
90 public:
91  typedef Opm::LiquidPhase<Scalar, Opm::DNAPL<Scalar> > type;
92 };
93 
94 // Set the material Law
95 SET_PROP(LensBaseProblem, MaterialLaw)
96 {
97 private:
98  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
99  enum { wettingPhaseIdx = FluidSystem::wettingPhaseIdx };
100  enum { nonWettingPhaseIdx = FluidSystem::nonWettingPhaseIdx };
101 
102  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
103  typedef Opm::TwoPhaseMaterialTraits<Scalar,
104  /*wettingPhaseIdx=*/FluidSystem::wettingPhaseIdx,
105  /*nonWettingPhaseIdx=*/FluidSystem::nonWettingPhaseIdx> Traits;
106 
107  // define the material law which is parameterized by effective
108  // saturations
109  typedef Opm::RegularizedVanGenuchten<Traits> EffectiveLaw;
110 
111 public:
112  // define the material law parameterized by absolute saturations
113  typedef Opm::EffToAbsLaw<EffectiveLaw> type;
114 };
115 
116 // Write the solutions of individual newton iterations?
117 SET_BOOL_PROP(LensBaseProblem, NewtonWriteConvergence, false);
118 
119 // Use forward differences instead of central differences
120 SET_INT_PROP(LensBaseProblem, NumericDifferenceMethod, +1);
121 
122 // Enable gravity
123 SET_BOOL_PROP(LensBaseProblem, EnableGravity, true);
124 
125 // define the properties specific for the lens problem
126 SET_SCALAR_PROP(LensBaseProblem, LensLowerLeftX, 1.0);
127 SET_SCALAR_PROP(LensBaseProblem, LensLowerLeftY, 2.0);
128 SET_SCALAR_PROP(LensBaseProblem, LensLowerLeftZ, 0.0);
129 SET_SCALAR_PROP(LensBaseProblem, LensUpperRightX, 4.0);
130 SET_SCALAR_PROP(LensBaseProblem, LensUpperRightY, 3.0);
131 SET_SCALAR_PROP(LensBaseProblem, LensUpperRightZ, 1.0);
132 
133 SET_SCALAR_PROP(LensBaseProblem, DomainSizeX, 6.0);
134 SET_SCALAR_PROP(LensBaseProblem, DomainSizeY, 4.0);
135 SET_SCALAR_PROP(LensBaseProblem, DomainSizeZ, 1.0);
136 
137 SET_INT_PROP(LensBaseProblem, CellsX, 48);
138 SET_INT_PROP(LensBaseProblem, CellsY, 32);
139 SET_INT_PROP(LensBaseProblem, CellsZ, 16);
140 
141 // The default for the end time of the simulation
142 SET_SCALAR_PROP(LensBaseProblem, EndTime, 30e3);
143 
144 // The default for the initial time step size of the simulation
145 SET_SCALAR_PROP(LensBaseProblem, InitialTimeStepSize, 250);
146 
147 // By default, include the intrinsic permeability tensor to the VTK output files
148 SET_BOOL_PROP(LensBaseProblem, VtkWriteIntrinsicPermeabilities, true);
149 
150 // enable the storage cache by default for this problem
151 SET_BOOL_PROP(LensBaseProblem, EnableStorageCache, true);
152 
153 // enable the cache for intensive quantities by default for this problem
154 SET_BOOL_PROP(LensBaseProblem, EnableIntensiveQuantityCache, true);
155 } // namespace Properties
156 
180 template <class TypeTag>
181 class LensProblem : public GET_PROP_TYPE(TypeTag, BaseProblem)
182 {
183  typedef typename GET_PROP_TYPE(TypeTag, BaseProblem) ParentType;
184 
185  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
186  typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
187  typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
188  typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
189  typedef typename GET_PROP_TYPE(TypeTag, WettingPhase) WettingPhase;
190  typedef typename GET_PROP_TYPE(TypeTag, NonwettingPhase) NonwettingPhase;
191  typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
192  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
193  typedef typename GET_PROP_TYPE(TypeTag, Model) Model;
194 
195  enum {
196  // number of phases
197  numPhases = FluidSystem::numPhases,
198 
199  // phase indices
200  wettingPhaseIdx = FluidSystem::wettingPhaseIdx,
201  nonWettingPhaseIdx = FluidSystem::nonWettingPhaseIdx,
202 
203  // equation indices
204  contiNEqIdx = Indices::conti0EqIdx + nonWettingPhaseIdx,
205 
206  // Grid and world dimension
207  dim = GridView::dimension,
208  dimWorld = GridView::dimensionworld
209  };
210 
211  typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
212  typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
213  typedef typename GET_PROP_TYPE(TypeTag, BoundaryRateVector) BoundaryRateVector;
214  typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
215  typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
216 
217  typedef typename GridView::ctype CoordScalar;
218  typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
219 
220  typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
221 
222 public:
226  LensProblem(Simulator& simulator)
227  : ParentType(simulator)
228  { }
229 
233  void finishInit()
234  {
235  ParentType::finishInit();
236 
237  eps_ = 3e-6;
238  FluidSystem::init();
239 
240  temperature_ = 273.15 + 20; // -> 20°C
241  lensLowerLeft_[0] = EWOMS_GET_PARAM(TypeTag, Scalar, LensLowerLeftX);
242  lensLowerLeft_[1] = EWOMS_GET_PARAM(TypeTag, Scalar, LensLowerLeftY);
243  lensUpperRight_[0] = EWOMS_GET_PARAM(TypeTag, Scalar, LensUpperRightX);
244  lensUpperRight_[1] = EWOMS_GET_PARAM(TypeTag, Scalar, LensUpperRightY);
245 
246  if (dimWorld == 3) {
247  lensLowerLeft_[2] = EWOMS_GET_PARAM(TypeTag, Scalar, LensLowerLeftZ);
248  lensUpperRight_[2] = EWOMS_GET_PARAM(TypeTag, Scalar, LensUpperRightZ);
249  }
250 
251  // residual saturations
252  lensMaterialParams_.setResidualSaturation(wettingPhaseIdx, 0.18);
253  lensMaterialParams_.setResidualSaturation(nonWettingPhaseIdx, 0.0);
254  outerMaterialParams_.setResidualSaturation(wettingPhaseIdx, 0.05);
255  outerMaterialParams_.setResidualSaturation(nonWettingPhaseIdx, 0.0);
256 
257  // parameters for the Van Genuchten law: alpha and n
258  lensMaterialParams_.setVgAlpha(0.00045);
259  lensMaterialParams_.setVgN(7.3);
260  outerMaterialParams_.setVgAlpha(0.0037);
261  outerMaterialParams_.setVgN(4.7);
262 
263  lensMaterialParams_.finalize();
264  outerMaterialParams_.finalize();
265 
266  lensK_ = this->toDimMatrix_(9.05e-12);
267  outerK_ = this->toDimMatrix_(4.6e-10);
268 
269  if (dimWorld == 3) {
270  this->gravity_ = 0;
271  this->gravity_[1] = -9.81;
272  }
273  }
274 
278  static void registerParameters()
279  {
280  ParentType::registerParameters();
281 
282  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensLowerLeftX,
283  "The x-coordinate of the lens' lower-left corner "
284  "[m].");
285  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensLowerLeftY,
286  "The y-coordinate of the lens' lower-left corner "
287  "[m].");
288  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensUpperRightX,
289  "The x-coordinate of the lens' upper-right corner "
290  "[m].");
291  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensUpperRightY,
292  "The y-coordinate of the lens' upper-right corner "
293  "[m].");
294 
295  if (dimWorld == 3) {
296  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensLowerLeftZ,
297  "The z-coordinate of the lens' lower-left "
298  "corner [m].");
299  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensUpperRightZ,
300  "The z-coordinate of the lens' upper-right "
301  "corner [m].");
302  }
303  }
304 
308 
313  template <class Context>
314  const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
315  unsigned timeIdx) const
316  {
317  const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
318 
319  if (isInLens_(globalPos))
320  return lensK_;
321  return outerK_;
322  }
323 
327  template <class Context>
328  Scalar porosity(const Context& context OPM_UNUSED,
329  unsigned spaceIdx OPM_UNUSED,
330  unsigned timeIdx OPM_UNUSED) const
331  { return 0.4; }
332 
336  template <class Context>
337  const MaterialLawParams& materialLawParams(const Context& context,
338  unsigned spaceIdx, unsigned timeIdx) const
339  {
340  const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
341 
342  if (isInLens_(globalPos))
343  return lensMaterialParams_;
344  return outerMaterialParams_;
345  }
346 
350  template <class Context>
351  Scalar temperature(const Context& context OPM_UNUSED,
352  unsigned spaceIdx OPM_UNUSED,
353  unsigned timeIdx OPM_UNUSED) const
354  { return temperature_; }
355 
357 
361 
366  std::string name() const
367  {
368  typedef typename GET_PROP_TYPE(TypeTag, LocalLinearizerSplice) LLS;
369 
370  bool useAutoDiff = std::is_same<LLS, TTAG(AutoDiffLocalLinearizer)>::value;
371 
372  std::ostringstream oss;
373  oss << "lens_" << Model::name()
374  << "_" << Model::discretizationName()
375  << "_" << (useAutoDiff?"ad":"fd");
376  return oss.str();
377  }
378 
383  { }
384 
389  { }
390 
394  void endTimeStep()
395  {
396 #ifndef NDEBUG
397  this->model().checkConservativeness();
398 
399  // Calculate storage terms
400  EqVector storage;
401  this->model().globalStorage(storage);
402 
403  // Write mass balance information for rank 0
404  if (this->gridView().comm().rank() == 0) {
405  std::cout << "Storage: " << storage << std::endl << std::flush;
406  }
407 #endif // NDEBUG
408  }
409 
411 
415 
420  template <class Context>
421  void boundary(BoundaryRateVector& values,
422  const Context& context,
423  unsigned spaceIdx,
424  unsigned timeIdx) const
425  {
426  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
427 
428  if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
429  // free flow boundary
430  Scalar densityW = WettingPhase::density(temperature_,
431  /*pressure=*/Scalar(1e5));
432 
433  Scalar T = temperature(context, spaceIdx, timeIdx);
434  Scalar pw, Sw;
435 
436  // set wetting phase pressure and saturation
437  if (onLeftBoundary_(pos)) {
438  Scalar height = this->boundingBoxMax()[1] - this->boundingBoxMin()[1];
439  Scalar depth = this->boundingBoxMax()[1] - pos[1];
440  Scalar alpha = (1 + 1.5 / height);
441 
442  // hydrostatic pressure scaled by alpha
443  pw = 1e5 - alpha * densityW * this->gravity()[1] * depth;
444  Sw = 1.0;
445  }
446  else {
447  Scalar depth = this->boundingBoxMax()[1] - pos[1];
448 
449  // hydrostatic pressure
450  pw = 1e5 - densityW * this->gravity()[1] * depth;
451  Sw = 1.0;
452  }
453 
454  // specify a full fluid state using pw and Sw
455  const MaterialLawParams& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
456 
457  Opm::ImmiscibleFluidState<Scalar, FluidSystem,
458  /*storeEnthalpy=*/false> fs;
459  fs.setSaturation(wettingPhaseIdx, Sw);
460  fs.setSaturation(nonWettingPhaseIdx, 1 - Sw);
461  fs.setTemperature(T);
462 
463  Scalar pC[numPhases];
464  MaterialLaw::capillaryPressures(pC, matParams, fs);
465  fs.setPressure(wettingPhaseIdx, pw);
466  fs.setPressure(nonWettingPhaseIdx, pw + pC[nonWettingPhaseIdx] - pC[wettingPhaseIdx]);
467 
468  // impose an freeflow boundary condition
469  values.setFreeFlow(context, spaceIdx, timeIdx, fs);
470  }
471  else if (onInlet_(pos)) {
472  RateVector massRate(0.0);
473  massRate = 0.0;
474  massRate[contiNEqIdx] = -0.04; // kg / (m^2 * s)
475 
476  // impose a forced flow boundary
477  values.setMassRate(massRate);
478  }
479  else {
480  // no flow boundary
481  values.setNoFlow();
482  }
483  }
484 
486 
490 
495  template <class Context>
496  void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx, unsigned timeIdx) const
497  {
498  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
499  Scalar depth = this->boundingBoxMax()[1] - pos[1];
500 
501  Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
502  fs.setPressure(wettingPhaseIdx, /*pressure=*/1e5);
503 
504  Scalar Sw = 1.0;
505  fs.setSaturation(wettingPhaseIdx, Sw);
506  fs.setSaturation(nonWettingPhaseIdx, 1 - Sw);
507 
508  fs.setTemperature(temperature_);
509 
510  typename FluidSystem::template ParameterCache<Scalar> paramCache;
511  paramCache.updatePhase(fs, wettingPhaseIdx);
512  Scalar densityW = FluidSystem::density(fs, paramCache, wettingPhaseIdx);
513 
514  // hydrostatic pressure (assuming incompressibility)
515  Scalar pw = 1e5 - densityW * this->gravity()[1] * depth;
516 
517  // calculate the capillary pressure
518  const MaterialLawParams& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
519  Scalar pC[numPhases];
520  MaterialLaw::capillaryPressures(pC, matParams, fs);
521 
522  // make a full fluid state
523  fs.setPressure(wettingPhaseIdx, pw);
524  fs.setPressure(nonWettingPhaseIdx, pw + (pC[wettingPhaseIdx] - pC[nonWettingPhaseIdx]));
525 
526  // assign the primary variables
527  values.assignNaive(fs);
528  }
529 
536  template <class Context>
537  void source(RateVector& rate,
538  const Context& context OPM_UNUSED,
539  unsigned spaceIdx OPM_UNUSED,
540  unsigned timeIdx OPM_UNUSED) const
541  { rate = Scalar(0.0); }
542 
544 
545 private:
546  bool isInLens_(const GlobalPosition& pos) const
547  {
548  for (unsigned i = 0; i < dim; ++i) {
549  if (pos[i] < lensLowerLeft_[i] - eps_ || pos[i] > lensUpperRight_[i]
550  + eps_)
551  return false;
552  }
553  return true;
554  }
555 
556  bool onLeftBoundary_(const GlobalPosition& pos) const
557  { return pos[0] < this->boundingBoxMin()[0] + eps_; }
558 
559  bool onRightBoundary_(const GlobalPosition& pos) const
560  { return pos[0] > this->boundingBoxMax()[0] - eps_; }
561 
562  bool onLowerBoundary_(const GlobalPosition& pos) const
563  { return pos[1] < this->boundingBoxMin()[1] + eps_; }
564 
565  bool onUpperBoundary_(const GlobalPosition& pos) const
566  { return pos[1] > this->boundingBoxMax()[1] - eps_; }
567 
568  bool onInlet_(const GlobalPosition& pos) const
569  {
570  Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0];
571  Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width;
572  return onUpperBoundary_(pos) && 0.5 < lambda && lambda < 2.0 / 3.0;
573  }
574 
575  GlobalPosition lensLowerLeft_;
576  GlobalPosition lensUpperRight_;
577 
578  DimMatrix lensK_;
579  DimMatrix outerK_;
580  MaterialLawParams lensMaterialParams_;
581  MaterialLawParams outerMaterialParams_;
582 
583  Scalar temperature_;
584  Scalar eps_;
585 };
586 
587 } // namespace Ewoms
588 
589 #endif
static void registerParameters()
Definition: lensproblem.hh:278
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
const DimMatrix & intrinsicPermeability(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: lensproblem.hh:314
Definition: baseauxiliarymodule.hh:37
Helper class for grid instantiation of the lens problem.
Definition: structuredgridmanager.hh:51
LensProblem(Simulator &simulator)
Definition: lensproblem.hh:226
void source(RateVector &rate, const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Evaluate the source term for all phases within a given sub-control-volume.
Definition: lensproblem.hh:537
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
Helper class for grid instantiation of the lens problem.
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
void boundary(BoundaryRateVector &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Evaluate the boundary conditions for a boundary segment.
Definition: lensproblem.hh:421
#define INHERITS_FROM(...)
Syntactic sugar for NEW_TYPE_TAG.
Definition: propertysystem.hh:230
#define SET_INT_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant integer value.
Definition: propertysystem.hh:345
void beginTimeStep()
Called by the simulator before each time integration.
Definition: lensproblem.hh:382
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:68
void beginIteration()
Called by the simulator before each Newton-Raphson iteration.
Definition: lensproblem.hh:388
void endTimeStep()
Called by the simulator after each time integration.
Definition: lensproblem.hh:394
std::string name() const
The problem name.
Definition: lensproblem.hh:366
Scalar temperature(const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Definition: lensproblem.hh:351
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
Calculates the local residual and its Jacobian for a single element of the grid.
void finishInit()
Called by the Ewoms::Simulator in order to initialize the problem.
Definition: lensproblem.hh:233
#define SET_PROP(EffTypeTagName, PropTagName)
Set a property for a specific type tag.
Definition: propertysystem.hh:297
Scalar porosity(const Context &context OPM_UNUSED, unsigned spaceIdx OPM_UNUSED, unsigned timeIdx OPM_UNUSED) const
Definition: lensproblem.hh:328
const MaterialLawParams & materialLawParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: lensproblem.hh:337
Defines the properties required for the immiscible multi-phase model.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
#define NEW_PROP_TAG(PTagName)
Define a property tag.
Definition: propertysystem.hh:247
#define SET_TYPE_PROP(EffTypeTagName, PropTagName,...)
Set a property which defines a type.
Definition: propertysystem.hh:377
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Evaluate the initial value for a control volume.
Definition: lensproblem.hh:496
Soil contamination problem where DNAPL infiltrates a fully water saturated medium.
Definition: lensproblem.hh:55
#define SET_SCALAR_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant scalar value.
Definition: propertysystem.hh:394