Computing Raviart Thomas function

Hi guys, can someone help me? I’m having a hard time computing the following set of equations:
Screenshot from 2023-07-02 12-13-48

This is the mesh I’m working with:

I need to compute the normal component of sigma on all egdes, I tried doing like:

func RT = RT0;
func PPme = P0edge;

fespace RTm(Th[K], RT);
RTm [sigmaHx, sigmaHy];
fespace PPm(Th[K], PPme);
PPm q;

matrix mCk1;

    varf chk1([sigmaHx, sigmaHy], [q])
      = intalledges(Th[K])((sigmaHx*N.x + sigmaHy*N.y)*q);
    mCk1 = chk1(RTm, PPm);
    cout << mCk1 << "\n";

but it does not seem to work on internal edges

for the RHS, I tried something like:

VhKaux AgraduHhX = kappa*dx(auxuHhK);
VhKaux AgraduHhY = kappa*dy(auxuHhK);
    varf lhK([vsolx, vsoly], [q])
      = intalledges(Th[K])((1-nTonEdge)*(mean(AgraduHhX)*N.x + mean(AgraduHhY)*N.y)*q);

    real[int] vlhK = lhK(0, PPm);
    // cout << vlhK << "\n";

    real[int] vbl(ndofPPm);
    for(int e=0; e<edgeLabel.n; e++){
      varf blk(muHk,q)
      = int1d(Th[K], edgeLabel[e])(edgeSignG(K,e)*lambdaH*q);

      real[int] vblaux = blk(0,PPm);
      vbl = vbl + vblaux;
      // cout << vblaux << "\n";
    // cout << vbl << "\n";

does anyone have any idea if this is right?