Add support for solving scalars using a different polynomial Order

[stgeke]

• [ ] Replace cds->gsh by cds->gsh[] + eliminate cds->gshT
• [ ] For cds kernels using p_Nq or p_subNq we have to use an array of kernel pointers!
• [ ] How do we ensure that the elliptic working are big enough for scalars (currently we’re using nrs->fieldOffset)
• [ ] How do we ensure that platform->mempool is big enough? (currently we’re using nrs->fieldOffset)
• [ ] Interpolate velocity onto cds->mesh[] if required
• [ ] Interpolate scalar to velocity N for output (we cannot have two different N in the same field file)