Store coefficient values at quadrature points for kernels

Aperta
#565 0 commenti 0 reazioni 0 assegnatari Vedi su GitHub

I maintainer di solito rispondono entro 1 giorno

Nessuno ha ancora preso questa issue.

Valutazione

Difficoltà
5/5
Tempo stimato
Più di una settimana
Idoneità per principianti
35/100
Tipo di issue
Funzionalità
Chiarezza
Abbastanza chiara
Stato di attività
Ferma
Stack tecnologico
cpp

Direzione di ricerca

Iniziate con il codice di assemblaggio del kernel illustrato nell’issue. Confrontate l’archiviazione corrente per punto con gli array di punti di quadratura proposti e individuate il percorso di pull-back e delle informazioni geometriche che rimane non vettorizzato. Il lavoro è completato quando i dati dei coefficienti e i dati geometrici sono memorizzati nei punti di quadratura, in modo che questo lavoro possa essere vettorizzato sui punti di quadratura senza modificare i risultati assemblati.

Scritto dal modello di indicizzazione a partire dal testo della issue.

Descrizione

Except for few modifications the code for assembling a vector looks like:

  T w1[3] = {0};
  T w0[1] = {0};
  for (int ip = 0; ip < Np; ip++)
  {
    // vectorized across
    for (int ic = 0; ic < Nc; ++ic)
    {
      w1[0] += w[4 + ic] * D0[ip][ic];
      w1[1] += w[4 + ic] * D1[ip][ic];
      w1[2] += w[4 + ic] * D2[ip][ic];
    }
    for (int ic = 0; ic < 4; ++ic)
      w0[0] += w[ic] * Phi[ip][ic];
    
    // Not vectorized
    T fw[3] = {0};
    fw[0] = (G[0][0] * w1[0] + G[0][1] * w1[1] + G[0][2] * w1[2]) * weights[ip];
    fw[1] = (G[1][0] * w1[0] + G[1][1] * w1[1] + G[1][2] * w1[2]) * weights[ip];
    fw[2] = (G[2][0] * w1[0] + G[2][1] * w1[1] + G[2][2] * w1[2]) * weights[ip];

    for (int i = 0; i < Nc; ++i)
      A[i] += fw[0 * Np + ip] * D0[ip][i] + fw[1 * Np + ip] * D1[ip][i] + fw[2 * Np + ip] * D2[ip][i];
  }

The application of pull-back and geometric information is not vectorized.
We should consider storing more data at quadrature points to enable vectorization across quadrature points.

    T w1[3 * Np] = {0};
    T w0[Np] = {0};
    for (int ip = 0; ip < Np; ip++)
    {
      for (int ic = 0; ic < Nc; ++ic)
      {
        w1[0 * Np + ip] += w[4 + ic] * D0[ip][ic];
        w1[1 * Np + ip] += w[4 + ic] * D1[ip][ic];
        w1[2 * Np + ip] += w[4 + ic] * D2[ip][ic];
      }
      for (int ic = 0; ic < 4; ++ic)
        w0[ip] += w[ic] * Phi[ip][ic];
    }

    T fw[3 * Np] = {0};
     // Vectorized across quadrature points
    for (int ip = 0; ip < Np; ip++)
    {
      fw[0 * Np + ip] = (G[0][0] * w1[ip] + G[0][1] * w1[Np + ip] + G[0][2] * w1[2 * Np + ip]) * weights[ip ];
      fw[1 * Np + ip] = (G[1][0] * w1[ip] + G[1][1] * w1[Np + ip] + G[1][2] * w1[2 * Np + ip]) * weights[ip ];
      fw[2 * Np + ip] = (G[2][0] * w1[ip] + G[2][1] * w1[Np + ip] + G[2][2] * w1[2 * Np + ip]) * weights[ip ];
    }
    for (int ip = 0; ip < Np; ip++)
    {
      for (int i = 0; i < Nc; ++i)
        A[i] += fw[0 * Np + ip] * D0[ip][i] + fw[1 * Np + ip] * D1[ip][i] + fw[2 * Np + ip] * D2[ip][i];
    }
Lingua principale
Python
Stelle
196
Fork
47
Merge medio
2g 13h
PR unite (30g)
24

Preparare l'ambiente

Questo progetto non fornisce container di sviluppo, Dockerfile né guida per i contributori, quindi l'ambiente è a tuo carico: parti dal suo README e consulta la nostra guida al primo contributo per i passaggi generali.

Come iniziare

  1. Leggi tutta la issue e poi la guida ai contributi del progetto.
  2. Commenta sulla issue per dire che te ne occupi tu — evita che due persone facciano lo stesso lavoro.
  3. Fai un fork del repository e lavora su un branch.
  4. Apri una pull request che faccia riferimento al numero della issue.

Altre issue di FEniCS/ffcx

Tutte le issue di FEniCS/ffcx

Issue simili

Altre issue su Python

Ricevi le nuove issue nella tua casella

Un breve riepilogo di issue GitHub adatte ai principianti.