if (_source_group.containsMyRank())
{
mcIdType nbrows = ToIdType(_coeffs.size());
-
// performing W.S
// W is the intersection matrix
// S is the source vector
if (_target_group.containsMyRank())
{
- mcIdType nbelems = field.getArray()->getNumberOfTuples() ;
- double* value = const_cast<double*> (field.getArray()->getPointer());
- for (mcIdType i=0; i<nbelems*nbcomp; i++)
- {
- value[i]=0.0;
- }
+ field.getArray()->fillWithZero();
}
//on source side : sending T=VT^(-1).(W.S)
//on target side :: receiving T and storing it in field
- _mapping.sendRecv(&target_value[0],field);
+ _mapping.sendRecv(target_value.data(),field,this->_presence_dft_value,this->_dft_value);
}
//T is the target vector
for (mcIdType irow = 0; irow < nbrows; irow++)
{
- for (mcIdType icol = _row_offsets[irow]; icol < _row_offsets[irow+1]; icol++)
- {
- int colid = _coeffs[irow][icol-_row_offsets[irow]].first;
- double value = _coeffs[irow][icol-_row_offsets[irow]].second;
- double deno = _deno_reverse_multiply[irow][icol-_row_offsets[irow]];
- for (std::size_t icomp=0; icomp<nbcomp; icomp++)
- {
- double coeff_row = source_value[colid*nbcomp+icomp];
- array[irow*nbcomp+icomp] += value*coeff_row/deno;
- }
- }
+ if( _row_offsets[irow+1] > _row_offsets[irow] )
+ {
+ for (mcIdType icol = _row_offsets[irow]; icol < _row_offsets[irow+1]; icol++)
+ {
+ int colid = _coeffs[irow][icol-_row_offsets[irow]].first;
+ double value = _coeffs[irow][icol-_row_offsets[irow]].second;
+ double deno = _deno_reverse_multiply[irow][icol-_row_offsets[irow]];
+ for (std::size_t icomp=0; icomp<nbcomp; icomp++)
+ {
+ double coeff_row = source_value[colid*nbcomp+icomp];
+ array[irow*nbcomp+icomp] += value*coeff_row/deno;
+ }
+ }
+ }
+ else
+ {
+ if( _presence_dft_value )
+ std::fill(array+irow*nbcomp,array+(irow+1)*nbcomp,this->_dft_value);
+ }
}
}
}
* The ids that were defined by addElementFromSource method
* are sent.
*/
- void MxN_Mapping::sendRecv(double* sendfield, MEDCouplingFieldDouble& field) const
+ void MxN_Mapping::sendRecv(double* sendfield, MEDCouplingFieldDouble& field, bool isDftValue, double dftValue) const
{
CommInterface comm_interface=_union_group->getCommInterface();
const MPIProcessorGroup* group = static_cast<const MPIProcessorGroup*>(_union_group);
sendbuf = new double[_sending_ids.size()*nbcomp];
if (_recv_ids.size()>0)
recvbuf = new double[_recv_ids.size()*nbcomp];
-
int* sendcounts = new int[_union_group->size()];
int* senddispls=new int[_union_group->size()];
int* recvcounts=new int[_union_group->size()];
//setting the received values in the field
DataArrayDouble *fieldArr=field.getArray();
- double* recvptr=recvbuf;
+ double* recvptr=recvbuf;
+ mcIdType nbTuples(fieldArr->getNumberOfTuples());
+ std::vector<bool> hitMachine(nbTuples,false);
for (int i=0; i< _recv_proc_offsets[_union_group->size()]; i++)
{
+ mcIdType recvId(_recv_ids[i]);
+ hitMachine[recvId] = true;
for (int icomp=0; icomp<nbcomp; icomp++)
{
- double temp = fieldArr->getIJ(_recv_ids[i],icomp);
- fieldArr->setIJ(_recv_ids[i],icomp,temp+*recvptr);
+ double temp = fieldArr->getIJ(recvId,icomp);
+ fieldArr->setIJ(recvId,icomp,temp+*recvptr);
recvptr++;
}
- }
+ }
+ if( isDftValue )
+ {
+ double *fieldPtr( fieldArr->getPointerSilent() );
+ for( mcIdType ituple = 0 ; ituple < nbTuples ; ++ituple )
+ {
+ if( ! hitMachine[ituple] )
+ std::fill( fieldPtr + ituple*nbcomp, fieldPtr + (ituple+1)*nbcomp, dftValue);
+ }
+ }
if (sendbuf!=0 && getAllToAllMethod()== Native)
delete[] sendbuf;
if (recvbuf !=0)
source_group.release()
MPI.COMM_WORLD.Barrier()
+ def test_InterpKernelDEC_default(self):
+ """
+ [EDF27375] : Put a default value when non intersecting case
+ """
+ size = MPI.COMM_WORLD.size
+ rank = MPI.COMM_WORLD.rank
+ if size != 4:
+ print("Should be run on 4 procs!")
+ return
+ nproc_source = 2
+ procs_source = list(range(nproc_source))
+ procs_target = list(range(size - nproc_source, size))
+
+ interface = CommInterface()
+ target_group = MPIProcessorGroup(interface, procs_target)
+ source_group = MPIProcessorGroup(interface, procs_source)
+ dec = InterpKernelDEC(source_group, target_group)
+ #
+ MPI.COMM_WORLD.Barrier()
+ if source_group.containsMyRank():
+ mesh = eval("Source_Proc_{}".format(rank))()
+ nb_local=mesh.getNumberOfCells()
+ field = MEDCouplingFieldDouble(ON_CELLS)
+ field.setNature(IntensiveMaximum)
+ field.setMesh(mesh)
+ arr = DataArrayDouble(nb_local) ; arr.iota() ; arr += rank
+ field.setArray(arr)
+ dec.attachLocalField(field)
+ dec.synchronizeWithDefaultValue(-2000.0)
+ dec.sendData()
+ # target -> source
+ dec.recvData()
+ if rank == 0:
+ self.assertTrue(field.getArray().isEqual(DataArrayDouble([0.6,0.6,-2000]),1e-12))
+ if rank == 1:
+ self.assertTrue(field.getArray().isEqual(DataArrayDouble([1.0,-2000]),1e-12))
+ else:
+ mesh = eval("Target_Proc_{}".format(rank))()
+ nb_local=mesh.getNumberOfCells()
+ field = MEDCouplingFieldDouble(ON_CELLS)
+ field.setNature(IntensiveMaximum)
+ field.setMesh(mesh)
+ arr = DataArrayDouble(nb_local) ; arr[:] = -20
+ field.setArray(arr)
+ dec.attachLocalField(field)
+ dec.synchronizeWithDefaultValue(-1000.0)
+ dec.recvData()
+ if rank == 2:
+ # matrix S0 / T2 = [[(0,S0,1),(1,S0,1.5)]
+ # IntensiveMaximum => [[(0,S0,1/2.5),(1,S0,1.5/2.5)]
+ #
+ self.assertTrue(field.getArray().isEqual(DataArrayDouble([0.6]),1e-12))
+ if rank == 3:
+ # matrix S1 / T3 = [[],[(0,S1,1.0)],[(0,S1,2.0)],[]]
+ # IntensiveMaximum => [[],[(0,S1,1.0/1.0)],[(0,S1,2.0/2.0)],[]]
+ self.assertTrue(field.getArray().isEqual(DataArrayDouble([-1000.0, 1.0, 1.0, -1000.0]),1e-8))
+ # target -> source
+ dec.sendData()
+
+ # Some clean up that still needs MPI communication, so to be done before MPI_Finalize()
+ dec.release()
+ target_group.release()
+ source_group.release()
+ MPI.COMM_WORLD.Barrier()
+
+def Source_Proc_0():
+ coo = DataArrayDouble([(0,2),(2,2),(4,2),(0,4),(2,4),(4,4),(0,6),(2,6)])
+ m = MEDCouplingUMesh("mesh",2) ; m.setCoords(coo) ; m.allocateCells()
+ m.insertNextCell(NORM_QUAD4,[0,1,4,3])
+ m.insertNextCell(NORM_QUAD4,[1,2,5,4])
+ m.insertNextCell(NORM_QUAD4,[3,4,7,6])
+ return m
+
+def Source_Proc_1():
+ coo = DataArrayDouble([(6,2),(8,2),(10,2),(6,4),(8,4),(10,4)])
+ m = MEDCouplingUMesh("mesh",2) ; m.setCoords(coo) ; m.allocateCells()
+ m.insertNextCell(NORM_QUAD4,[0,1,4,3])
+ m.insertNextCell(NORM_QUAD4,[1,2,5,4])
+ return m
+
+def Target_Proc_2():
+ coo = DataArrayDouble([(1,0),(3.5,0),(1,3),(3.5,3)])
+ m = MEDCouplingUMesh("mesh",2) ; m.setCoords(coo) ; m.allocateCells()
+ m.insertNextCell(NORM_QUAD4,[0,1,3,2])
+ return m
+
+def Target_Proc_3():
+ coo = DataArrayDouble([(6,0),(7,0),(8,0),(9,0),(10,0),
+ (6,1),(7,1),(9,1),(10,1),
+ (7,3),(8,3),
+ (6,4),(7,4)])
+ m = MEDCouplingUMesh("mesh",2) ; m.setCoords(coo) ; m.allocateCells()
+ m.insertNextCell(NORM_QUAD4,[0,1,6,5])
+ m.insertNextCell(NORM_QUAD4,[1,2,10,9])
+ m.insertNextCell(NORM_QUAD4,[5,6,12,11])
+ m.insertNextCell(NORM_QUAD4,[3,4,8,7])
+ return m
+
if __name__ == "__main__":
unittest.main()
MPI.Finalize()
