Yes, that means that a significant part of MatShell logic would become default implementations. But I would now do that only with MatDiagonalScale as a proof-of-concept - then we can progressively do the same with the others. At some point it could happen even with axpy - I think it would be cool although more a complex task as you say. Vaclav
1. 12. 2018 v 21:31, Smith, Barry F. <[email protected]>:
So basically hoist most of the data structure
typedef struct { struct _MatShellOps ops[1];
PetscScalar vscale,vshift; Vec dshift; Vec left,right; Vec left_work,right_work; Vec left_add_work,right_add_work; Mat axpy; PetscScalar axpy_vscale; PetscBool managescalingshifts; /* The user will manage the scaling and shifts for the MATSHELL, not the default */ void *ctx; } Mat_Shell;
into Mat? MatShell handles the axpy (which is nasty), would you handle that as well, if not then you still have code duplication between the MatShell implementations and the Mat___Basic implementations. If it wasn't for the axpy I would say it is worth trying but including correctly the axpy (I am not sure if it is always handled correctly for MatShell) requires complex code.
Barry
On Dec 1, 2018, at 1:49 PM, Hapla Vaclav <[email protected]> wrote:
MatDiagonalScale_Shell, MatDiagonalScale_Normal, MatDiagonalScaleHermitian_Normal and perhaps some others are essentially the same.
What about converting them into one MatDiagonalScale_Basic which would be the default implementation if no type-specific implementation is provided?
I'm asking because I wanted to add the same to MATTRANSPOSE. But with my proposal, this redundancy would be removed and additionally MatDiagonalScale() would work for any current or future MatType while it doesn't prevent providing an optimised version for explicit matrix types.
Note the same could be done for MatScale, MatShift, MatDiagonalSet and maybe some others.
I also think such approach would promote using such high-level API functions in higher PETSc layers and user codes because one wouldn't have to be afraid those methods might not be implemented in input matrices.
Thanks for opinions
Vaclav