I just found EModelDefaultPrescribedEquationNumbering() which I hadn't seen before. Perfect! There are just some issues;
1. Construction of the element stiffness has to be repeated.
2. I could only get the K_bb and K_ii parts, as K_bi and K_ib requires a mix.
So, I will need to add another assembling routine for mixed UnknownNumberingScheme's.
The sparse matrices already support mixed assembling, so it's just one function to the engineering model.
I see that there is another assemble routine for mixed equation id. But i can't say I understand it's use.

Also, some small change needs to be done for the routines "buildInternalStructure" as it enforces square matrices.

Hi Borek.

The Kbb part is necessary to calculate the macroscopic tangent, even if it's not needed to solve the microscopic problem itself.
I realized also, that even if i had the whole matrix, splitting up the sparse matrices seemed like a very bad idea.

I can't see any reasonable way to compute this
C'*(K_bb - K_bi*K_ii\K_ib)*C
without obtaining these matrices induvidually. (C is a matrix i can construct knowing the prescribed dofs).

I didn't have to create any custom numbering, In fact, OOFEM already supports what i want to achieve, except for a few things;

* Building the internal structure for sparse matrices, it assumes to much about the size.
I implemented another version which respect the numberingscheme even for prescribed values, i.e;
buildInternalStructure(EngngModel *eModel, int di, EquationID ut, const UnknownNumberingScheme &r_s, const UnknownNumberingScheme &c_s)

* EngngModel->assemble didn't support mixed numberingschemes. Adding this was trivial.

* Petsc sparse matrices didn't support times(FloatArray &x,FloatArray &y);. But this was quite easy, so i implemented it already. I also added timesT (transpose) since I needed it as well, but I see no harm in that.

All other routines already natively supported

This seems to fit well into the design of OOFEM, and I don't think anything I've implemented would count as a hack, so I'm quite pleased with the results.

Since I need to do a bit of matrix operations (especially on sparse structures) I'm planning on filling up the gaps in the current implementations. Also for dense matrices, FloatMatrix and FloatArray. I'd like to have some basic BLAS operations.

The current code in buildInternalStructure assumes that the matrix is always square, and the size is the number of free equations, regardless of what numberingsscheme is submitted, so something definitely needs to be changed.

I dont think I could mix things up with another numbering scheme.
When it comes down to it, regardless of numbering scheme, EngngModel::assemble is just calling

        if ( mat.isNotEmpty() ) {
            if ( answer->assemble(loc, mat) == 0 ) {
                _error("assemble: sparse matrix assemble error");
            }
        }

with a single IntArray which is used for both columns and rows.

I'm suggesting that one adds the function

void EngngModel :: assemble(SparseMtrx *answer, TimeStep *tStep, EquationID ut,
                            CharType type, const UnknownNumberingScheme &rs, const UnknownNumberingScheme &cs,
                            Domain *domain)
// Same as assemble, but with different numbering for rows and columns
{
    int ielem;
    IntArray r_loc, c_loc;
    FloatMatrix mat;
    Element *element;

    if ( answer == NULL ) {
        _error("assemble: NULL pointer encountered.");
    }

    this->timer.resumeTimer(EngngModelTimer :: EMTT_NetComputationalStepTimer);
    int nelem = domain->giveNumberOfElements();
    for ( ielem = 1; ielem <= nelem; ielem++ ) {
        element = domain->giveElement(ielem);
#ifdef __PARALLEL_MODE
        // skip remote elements (these are used as mirrors of remote elements on other domains
        // when nonlocal constitutive models are used. They introduction is necessary to
        // allow local averaging on domains without fine grain communication between domains).
        if ( element->giveParallelMode() == Element_remote ) {
            continue;
        }

#endif

        this->giveElementCharacteristicMatrix(mat, ielem, type, tStep, domain);
        if ( mat.isNotEmpty() ) {
///////////////// I added this part, mixed numbering   ///////////////////////////
            element->giveLocationArray(r_loc, ut, rs);
            element->giveLocationArray(c_loc, ut, cs);
            if ( answer->assemble(r_loc, c_loc, mat) == 0 ) {
                _error("assemble: sparse matrix assemble error");
            }
        }
    }

    this->timer.pauseTimer(EngngModelTimer :: EMTT_NetComputationalStepTimer);

    answer->assembleBegin();
    answer->assembleEnd();
}

which is very similar the the routine for mixed EquationID, but for numberingschemes instead. I dont see how you could create a mixed matrix like this without adding this assemble routine.

I did this, and now I can assemble the mixed matrices by doing

    this->assemble( this->K_cf, tStep, EID_MomentumBalance_ConservationEquation,
            StiffnessMatrix, EModelDefaultPrescribedEquationNumbering(),
            EModelDefaultEquationNumbering(), this->giveDomain(1));

Nice to hear. I have also been filling in some gaps in the FEI classes. I will be starting to send you some patches by next week.