40 lduMatrix::preconditioner::
41 addasymMatrixConstructorToTable<distributedDILUPreconditioner>
63 const auto& lduAddr = matrix.
lduAddr();
73 for (
const label interfacei : selectedInterfaces)
75 interfaces[interfacei].initInterfaceMatrixUpdate
82 coupleCoeffs[interfacei],
91 for (
const label interfacei : selectedInterfaces)
93 interfaces[interfacei].updateInterfaceMatrix
100 coupleCoeffs[interfacei],
117 const auto& interfaces = solver_.interfaces();
119 if (selectedInterfaces.size())
122 FieldField<Field, scalar> one(interfaces.size());
123 FieldField<Field, solveScalar> old(interfaces.size());
124 for (
const label inti : selectedInterfaces)
126 const auto& intf = interfaces[inti].interface();
127 const auto& fc = intf.faceCells();
131 updateMatrixInterfaces
141 if (!colourBufs_.set(colouri))
146 new FieldField<Field, solveScalar>
152 auto& colourBuf = colourBufs_[colouri];
153 colourBuf.setSize(interfaces.size());
154 for (
const label inti : selectedInterfaces)
156 const auto& intf = interfaces[inti].interface();
157 const auto& fc = intf.faceCells();
158 if (!colourBuf.set(inti))
160 colourBuf.set(inti,
new Field<solveScalar>(fc.size()));
162 auto& cb = colourBuf[inti];
164 auto& oldValues = old[inti];
168 const label cell = fc[face];
170 cb[face] =
psi[cell]-oldValues[face];
172 std::swap(
psi[cell], oldValues[face]);
182 DynamicList<UPstream::Request>& requests
185 const auto& interfaces = solver_.interfaces();
186 const auto& interfaceBouCoeffs = solver_.interfaceBouCoeffs();
189 for (
const label inti : selectedInterfaces)
191 const auto& intf = interfaces[inti].interface();
192 const auto* ppp = isA<processorLduInterface>(intf);
194 auto& recvBuf = recvBufs_[inti];
197 requests.push_back(UPstream::Request());
202 recvBuf.data_bytes(),
203 recvBuf.size_bytes(),
215 DynamicList<UPstream::Request>& requests
218 const auto& interfaces = solver_.interfaces();
221 for (
const label inti : selectedInterfaces)
223 const auto& intf = interfaces[inti].interface();
224 const auto* ppp = isA<processorLduInterface>(intf);
225 const auto& faceCells = intf.faceCells();
227 auto& sendBuf = sendBufs_[inti];
229 sendBuf.resize_nocopy(faceCells.size());
232 sendBuf[face] = psiInternal[faceCells[face]];
235 requests.push_back(UPstream::Request());
240 sendBuf.cdata_bytes(),
241 sendBuf.size_bytes(),
251 DynamicList<UPstream::Request>& requests,
277 const auto& interfaces = solver_.interfaces();
278 const auto& interfaceBouCoeffs = solver_.interfaceBouCoeffs();
279 const auto& interfaceIntCoeffs = solver_.interfaceIntCoeffs();
281 const auto& intf = interfaces[inti].interface();
283 const auto&
faceCells = intf.faceCells();
284 const auto& bc = interfaceBouCoeffs[inti];
285 const auto& ic = interfaceIntCoeffs[inti];
304 const auto& matrix = solver_.matrix();
305 const auto& lduAddr = matrix.lduAddr();
307 const label*
const __restrict__ uPtr = lduAddr.upperAddr().begin();
308 const label*
const __restrict__ lPtr = lduAddr.lowerAddr().begin();
309 const scalar*
const __restrict__ upperPtr = matrix.upper().begin();
310 const scalar*
const __restrict__ lowerPtr = matrix.lower().begin();
312 const label nFaces = matrix.upper().size();
315 const auto& cellColour = *cellColourPtr_;
316 for (label face=0; face<nFaces; face++)
318 const label cell = lPtr[face];
319 if (cellColour[cell] == colouri)
321 rD[uPtr[face]] -= upperPtr[face]*lowerPtr[face]/rD[cell];
327 for (label face=0; face<nFaces; face++)
339 const Field<solveScalar>& recvBuf
342 const auto& interfaces = solver_.interfaces();
343 const auto& interfaceBouCoeffs = solver_.interfaceBouCoeffs();
345 const auto& intf = interfaces[inti].interface();
346 const auto& faceCells = intf.faceCells();
347 const auto& bc = interfaceBouCoeffs[inti];
348 const solveScalar* __restrict__ rDPtr = rD_.begin();
349 solveScalar* __restrict__ wAPtr = wA.begin();
366 const auto& matrix = solver_.matrix();
367 const auto& lduAddr = matrix.lduAddr();
369 solveScalar* __restrict__ wAPtr = wA.begin();
370 const solveScalar* __restrict__ rDPtr = rD_.begin();
372 const label*
const __restrict__ uPtr = lduAddr.upperAddr().begin();
373 const label*
const __restrict__ lPtr = lduAddr.lowerAddr().begin();
374 const label*
const __restrict__ losortPtr = lduAddr.losortAddr().begin();
375 const scalar*
const __restrict__ lowerPtr = matrix.lower().begin();
377 const label nFaces = matrix.upper().size();
380 const auto& cellColour = *cellColourPtr_;
381 for (label face=0; face<nFaces; face++)
383 const label sface = losortPtr[face];
384 const label cell = lPtr[sface];
385 if (cellColour[cell] == colouri)
387 wAPtr[uPtr[sface]] -=
388 rDPtr[uPtr[sface]]*lowerPtr[sface]*wAPtr[cell];
394 for (label face=0; face<nFaces; face++)
396 const label sface = losortPtr[face];
397 wAPtr[uPtr[sface]] -=
398 rDPtr[uPtr[sface]]*lowerPtr[sface]*wAPtr[lPtr[sface]];
410 const auto& matrix = solver_.matrix();
411 const auto& lduAddr = matrix.lduAddr();
413 solveScalar* __restrict__ wAPtr = wA.begin();
414 const solveScalar* __restrict__ rDPtr = rD_.begin();
416 const label*
const __restrict__ uPtr = lduAddr.upperAddr().begin();
417 const label*
const __restrict__ lPtr = lduAddr.lowerAddr().begin();
418 const scalar*
const __restrict__ upperPtr = matrix.upper().begin();
420 const label nFacesM1 = matrix.upper().size() - 1;
424 const auto& cellColour = *cellColourPtr_;
425 for (label face=nFacesM1; face>=0; face--)
427 const label cell = uPtr[face];
428 if (cellColour[cell] == colouri)
432 rDPtr[lPtr[face]]*upperPtr[face]*wAPtr[cell];
438 for (label face=nFacesM1; face>=0; face--)
452 const auto& matrix = solver_.matrix();
455 wait(lowerRecvRequests_);
458 receive(lowerNbrs_, lowerRecvRequests_);
462 rD.resize_nocopy(
diag.size());
463 std::copy(
diag.begin(),
diag.end(), rD.begin());
469 wait(lowerRecvRequests_);
471 for (
const label inti : lowerNbrs_)
473 addInterfaceDiag(rD, inti, recvBufs_[inti]);
487 for (label colouri = 0; colouri < nColours_; colouri++)
491 for (
const label inti : lowerGlobalRecv_[colouri])
493 addInterfaceDiag(rD, inti, colourBufs_[colouri][inti]);
497 forwardInternalDiag(rD, colouri);
502 sendGlobal(higherGlobalSend_[colouri], rD, higherColour_[colouri]);
507 wait(higherSendRequests_);
510 send(higherNbrs_, rD, higherSendRequests_);
513 const label nCells = rD.size();
515 for (label cell=0; cell<nCells; cell++)
517 rD[cell] = 1.0/rD[cell];
526 const lduMatrix::solver& sol,
527 const dictionary&
dict 530 lduMatrix::preconditioner(sol),
531 coupled_(
dict.getOrDefault<bool>(
"coupled", true, keyType::LITERAL)),
532 rD_(sol.matrix().
diag().size())
534 const lduMesh&
mesh = sol.matrix().mesh();
535 const auto& interfaces = sol.interfaces();
536 const auto& interfaceBouCoeffs = sol.interfaceBouCoeffs();
543 forAll(interfaceBouCoeffs, inti)
545 if (interfaces.set(inti))
547 const auto& bc = interfaceBouCoeffs[inti];
562 const label nProcColours =
567 Info<< typeName <<
" : calculated processor colours:" 568 << nProcColours <<
endl;
577 bool haveGlobalCoupled =
false;
578 bool haveDistributedAMI =
false;
581 if (interfaces.set(inti))
583 const auto& intf = interfaces[inti].interface();
584 const auto* ppp = isA<processorLduInterface>(intf);
587 const label nbrColour = procColours[ppp->neighbProcNo()];
589 if (nbrColour < myColour)
593 else if (nbrColour > myColour)
600 <<
"weird processorLduInterface" 601 <<
" from " << ppp->myProcNo()
602 <<
" to " << ppp->neighbProcNo()
608 haveGlobalCoupled =
true;
610 const auto* AMIpp = isA<cyclicAMILduInterface>(intf);
623 Info<< typeName <<
" : interface " << inti
624 <<
" of type " << intf.type()
625 <<
" is distributed:" << haveDistributedAMI <<
endl;
679 if (interfaces.set(inti))
681 const auto& intf = interfaces[inti].interface();
684 const auto* AMIpp = isA<cyclicAMILduInterface>(intf);
687 nbrInti = AMIpp->neighbPatchID();
689 const auto* cycpp = isA<cyclicLduInterface>(intf);
692 nbrInti = cycpp->neighbPatchID();
697 const label colouri = patchToColour[inti];
698 const label nbrColouri = patchToColour[nbrInti];
702 (haveDistributedAMI && (inti < nbrInti))
703 || (!haveDistributedAMI && (colouri < nbrColouri))
733 if (lowerRecv.size())
738 <<
"Lower global interfaces for colour:" << colouri
739 <<
" receiving from:" <<
flatOutput(lowerRecv)
741 <<
" with colour:" << lowerCol <<
endl;
747 if (higherRecv.size())
752 <<
"Higher global interfaces for colour:" << colouri
753 <<
" receiving from:" <<
flatOutput(higherRecv)
755 <<
" with colour:" << higherCol <<
endl;
772 wait(lowerSendRequests_);
773 wait(lowerRecvRequests_);
774 wait(higherSendRequests_);
775 wait(higherRecvRequests_);
794 solveScalar* __restrict__ wAPtr = wA.
begin();
795 const solveScalar* __restrict__ rAPtr = rA.
begin();
796 const solveScalar* __restrict__ rDPtr = rD_.begin();
798 const label nCells = wA.
size();
805 wait(lowerRecvRequests_);
808 receive(lowerNbrs_, lowerRecvRequests_);
819 wait(lowerRecvRequests_);
821 for (
const label inti : lowerNbrs_)
823 addInterface(wA, inti, recvBufs_[inti]);
826 for (label colouri = 0; colouri < nColours_; colouri++)
831 for (
const label inti : lowerGlobalRecv_[colouri])
833 addInterface(wA, inti, colourBufs_[colouri][inti]);
837 forwardInternal(wA, colouri);
845 higherGlobalSend_[colouri],
847 higherColour_[colouri]
854 wait(higherSendRequests_);
857 send(higherNbrs_, wA, higherSendRequests_);
864 wait(higherRecvRequests_);
867 receive(higherNbrs_, higherRecvRequests_);
871 wait(higherRecvRequests_);
873 for (
const label inti : higherNbrs_)
875 addInterface(wA, inti, recvBufs_[inti]);
878 for (label colouri = nColours_-1; colouri >= 0; colouri--)
883 for (
const label inti : higherGlobalRecv_[colouri])
885 addInterface(wA, inti, colourBufs_[colouri][inti]);
889 backwardInternal(wA, colouri);
897 lowerGlobalSend_[colouri],
899 lowerColour_[colouri]
906 wait(lowerSendRequests_);
909 send(lowerNbrs_, wA, lowerSendRequests_);
922 wait(lowerSendRequests_);
923 wait(lowerRecvRequests_);
924 wait(higherSendRequests_);
925 wait(higherRecvRequests_);
virtual void addInterfaceDiag(solveScalarField &rD, const label inti, const Field< solveScalar > &recvBuf) const
Update diagonal for interface.
FieldField< Field, solveScalar > sendBufs_
Buffers for sending and receiving data.
void size(const label n)
Older name for setAddressableSize.
static label cellColour(const lduMesh &mesh, labelList &cellColour, labelList &patchToColour)
Calculate (locally) per cell colour according to walking from global patches. Returns number of colou...
virtual void precondition(solveScalarField &wA, const solveScalarField &rA, const direction cmpt=0) const
Return wA the preconditioned form of residual rA.
Ostream & indent(Ostream &os)
Indent stream.
void send(const labelList &selectedInterfaces, const solveScalarField &psiInternal, DynamicList< UPstream::Request > &requests) const
Start sending sendBufs_.
static int incrMsgType(int val=1) noexcept
Increment the message tag for standard messages.
A face is a list of labels corresponding to mesh vertices.
Field< solveScalar > solveScalarField
distributedDILUPreconditioner(const lduMatrix::solver &, const dictionary &solverControlsUnused)
Construct from matrix components and preconditioner solver controls.
void append(const T &val)
Append an element at the end of the list.
const solver & solver_
Reference to the base-solver this preconditioner is used with.
static label nRequests() noexcept
Number of outstanding requests (on the internal list of requests)
autoPtr< labelList > cellColourPtr_
Local (cell) colouring from global interfaces.
const lduMatrix & matrix() const noexcept
virtual void forwardInternalDiag(solveScalarField &rD, const label colouri) const
Update diagonal for all faces of a certain colour.
solveScalarField rD_
The reciprocal preconditioned diagonal.
List< DynamicList< label > > lowerGlobalSend_
Interfaces to non-processor lower coupled interfaces.
virtual void setFinished(const solverPerformance &perf) const
Signal end of solver.
static const lduMesh * meshPtr_
Processor interface buffers and colouring.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual ~distributedDILUPreconditioner()
Destructor.
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
List< DynamicList< label > > higherGlobalSend_
Interfaces to non-processor higher coupled interfaces.
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
static int & msgType() noexcept
Message tag of standard messages.
static int myProcNo(const label communicator=worldComm)
Rank of this process in the communicator (starting from masterNo()). Can be negative if the process i...
Smooth ATC in cells next to a set of patches supplied by type.
FieldField< Field, solveScalar > recvBufs_
static void waitRequests()
Wait for all requests to finish.
A field of fields is a PtrList of fields with reference counting.
#define forAll(list, i)
Loop across all elements in list.
List< label > higherColour_
Corresponding destination colour (for higherGlobal)
static autoPtr< labelList > procColoursPtr_
Previous processor colours.
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
static label nProcs(const label communicator=worldComm)
Number of ranks in parallel run (for given communicator). It is 1 for serial run. ...
void setSize(const label n)
Alias for resize()
virtual void forwardInternal(solveScalarField &wA, const label colouri) const
Update preconditioned variable walking forward on internal faces.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual void addInterface(solveScalarField &wA, const label inti, const Field< solveScalar > &recvBuf) const
Update preconditioned variable from interface.
const lduInterfaceFieldPtrsList & interfaces() const noexcept
void append(const T &val)
Copy append an element to the end of this list.
Version of DILUpreconditioner that uses preconditioning across processor (and coupled) boundaries...
virtual void backwardInternal(solveScalarField &wA, const label colouri) const
Update preconditioned variable walking backward on internal faces.
iterator begin() noexcept
Return an iterator to begin traversing the UList.
virtual void calcReciprocalD(solveScalarField &rD) const
Calculate reciprocal of diagonal.
void receive(const labelList &selectedInterfaces, DynamicList< UPstream::Request > &requests) const
Start receiving in recvBufs_.
const bool coupled_
Precondition across global coupled bc.
void updateMatrixInterfaces(const bool add, const FieldField< Field, scalar > &coupleCoeffs, const labelList &selectedInterfaces, const solveScalarField &psiif, solveScalarField &result, const direction cmpt) const
Variant of lduMatrix::updateMatrixInterfaces on selected interfaces.
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
DynamicList< label > lowerNbrs_
Interfaces to lower coloured processors.
int debug
Static debugging option.
PtrList< FieldField< Field, solveScalar > > colourBufs_
Global interfaces. Per colour the interfaces that (might) influence it.
defineTypeNameAndDebug(combustionModel, 0)
static label colour(const lduMesh &mesh, labelList &procColour)
Calculate processor colouring from processor connectivity. Sets colour per processor such that no nei...
List< DynamicList< label > > lowerGlobalRecv_
Interfaces to non-processor lower coupled interfaces.
static void cancelRequests(UList< UPstream::Request > &requests)
Non-blocking comms: cancel and free outstanding requests. Corresponds to MPI_Cancel() + MPI_Request_f...
void sendGlobal(const labelList &selectedInterfaces, solveScalarField &psi, const label colouri) const
Send (and store in colourBufs_[colouri]) the effect of.
DynamicList< label > higherNbrs_
Interfaces to higher coloured processors.
label nColours_
Number of colours (in case of multiple disconnected regions.
#define WarningInFunction
Report a warning using Foam::Warning.
const lduAddressing & lduAddr() const
Return the LDU addressing.
A cell is defined as a list of faces with extra functionality.
static bool write(const UPstream::commsTypes commsType, const int toProcNo, const char *buf, const std::streamsize bufSize, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm, UPstream::Request *req=nullptr, const UPstream::sendModes sendMode=UPstream::sendModes::normal)
Write buffer contents to given processor.
#define DebugPout
Report an information message using Foam::Pout.
"nonBlocking" (immediate) : (MPI_Isend, MPI_Irecv)
messageStream Info
Information stream (stdout output on master, null elsewhere)
void wait(DynamicList< UPstream::Request > &requests, const bool cancel=false) const
Wait for requests or cancel/free requests.
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
const volScalarField & psi
List< label > labelList
A List of labels.
List< label > lowerColour_
Corresponding destination colour (for lowerGlobal)
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Ostream & incrIndent(Ostream &os)
Increment the indent level.
List< DynamicList< label > > higherGlobalRecv_
Interfaces to non-processor higher coupled interfaces.
lduMatrix::preconditioner::addasymMatrixConstructorToTable< distributedDILUPreconditioner > adddistributedDILUPreconditionerAsymMatrixConstructorToTable_
static std::streamsize read(const UPstream::commsTypes commsType, const int fromProcNo, char *buf, const std::streamsize bufSize, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm, UPstream::Request *req=nullptr)
Read buffer contents from given processor.
FlatOutput::OutputAdaptor< Container, Delimiters > flatOutput(const Container &obj, Delimiters delim)
Global flatOutput() function with specified output delimiters.
static constexpr const zero Zero
Global zero (0)