Champ (Field/SubField)

pyrucast a deux familles de champs : le champ aux nœuds (NodeField), qui donne une valeur à chaque nœud, et le champ aux points de Gauss (ElementField), qui en donne une à chaque point d’intégration de chaque maille. Les deux se lisent et se calculent de la même façon, et cette page décrit ce qu’elles ont en commun ; les deux pages suivantes donnent ce qui est propre à chacune.

En Python, il n’y a pas de classe Field : tout ce qui suit s’écrit à l’identique sur NodeField, sur ElementField et sur leurs zones (SubNodeField, SubElementField). Côté Rust, ce socle commun est écrit une seule fois, dans deux traits (voir Côté Rust).

Zones, supports et composantes

Un champ est un agrégat de zones. Chaque zone :

• est posée sur un support : un nuage de nœuds (une zone de maillage POI1) pour un champ aux nœuds, un sous-espace d’espace éléments finis pour un champ aux points de Gauss. zone.support() le rend, et champ.support() rend celui de chaque zone, dans l’ordre ;
• porte des composantes nommées ("T", "u_x", "E", "sigma_xx"…) : au moins une, chaque nom une seule fois ;
• range une valeur par composante et par ligne de son support : un nœud, ou un couple (maille, point de Gauss).

Les règles de cohérence

Tout champ respecte deux règles, que chaque opération garde :

1. Une composante par support. Sur un même support, une composante est portée par une seule zone. Plusieurs zones peuvent partager un support si leurs composantes diffèrent : les matériaux de deux physiques sur un même maillage, par exemple.
2. Une valeur par nœud. Un nœud que plusieurs zones portent y a la même valeur, pour chaque composante qu’elles ont en commun. Deux zones de maillage qui se touchent partagent leurs nœuds de frontière : le champ y a une valeur, rangée dans chacune. La règle ne concerne que le champ aux nœuds, car deux sous-espaces n’ont aucun point de Gauss en commun.

La première se vérifie à l’entrée. NodeField(maillage, …) refuse un maillage qui liste deux fois la même zone, et ElementField(espace, …) un espace qui liste deux fois le même sous-espace : l’union maillage | maillage garde chaque zone une fois. add_sub refuse une zone qui porterait une composante déjà portée sur son support.

La seconde se vérifie à l’union, et en sortie d’une opération entre deux champs quand c’est nécessaire (plus bas). Seule l’écriture dans une zone peut encore la rompre : champ[i].set_value(…) ne touche pas les autres zones qui portent le même nœud. Sur un champ aux nœuds, champ.check() la vérifie à la demande.

Lire un champ

Les zones d’un même champ n’ont pas forcément les mêmes composantes. L’exemple coupe un carré en deux zones de maillage, gauche et droite, qui partagent la diagonale : la température T est connue sur les deux, le flux q sur celle de droite seulement.

# T sur les deux zones, q sur celle de droite seulement.
t = pyrucast.NodeField(gauche, ["T"])
for noeud, valeur in zip(n[:3], (0.0, 10.0, 20.0)):
    t[0].set_value(noeud, "T", valeur)
tq = pyrucast.NodeField(droite, ["T", "q"])
for noeud, valeur in zip((n[1], n[3], n[2]), (10.0, 30.0, 20.0)):
    tq[0].set_value(noeud, "T", valeur)
tq[0].add_to_component("q", 5.0)
champ = t | tq

print(champ.components())  # ['T', 'q'] : les composantes de toutes les zones
print(champ.min("T"), champ.max("T"))  # 0.0 30.0 : cherchés dans les deux zones
print(champ.min())  # 0.0 : toutes les valeurs, toutes composantes confondues
print(champ.sum("q"))  # 15.0 : les trois nœuds de droite, à 5
print(len(champ["q"]), champ["q"].components())  # 1 ['q'] : la zone sans q disparaît

Une composante qu’aucune zone ne porte lève une erreur.

champ["c"] laisse de côté une zone qui ne porte aucune des composantes demandées. Une zone qui n’en porte pas d’autres est partagée avec le champ d’origine, sans copie : écrire dans l’une modifie l’autre. Une zone dont il faut retirer des composantes est recopiée. rename_component partage de même les zones qui ne portent pas a.

sum(c) ajoute une valeur par nœud : un nœud que deux zones partagent n’a qu’une valeur (deuxième règle), et compte une fois. Ce n’est pas une intégrale, qui pèse chaque valeur par la taille des mailles : c’est pyrucast.measure.integral. Deux autres réductions s’appliquent aux deux familles, et comptent de même un nœud une fois : pyrucast.measure.xtx(champ), la somme des carrés des valeurs, et pyrucast.measure.xty(a, b), la somme des produits des valeurs que a et b portent sur les mêmes supports et sous les mêmes noms (voir Opérateurs sur les champs).

Calculer avec des champs

Par un scalaire

champ + s, champ - s, champ * s, champ / s et champ ** s appliquent l’opération à chaque valeur, de chaque composante, de chaque zone, et rendent un nouveau champ ; le scalaire se place d’un côté ou de l’autre. -champ change le signe de chaque valeur. Il en va de même sur une zone.

plus_chaud = champ + 5.0  # un nouveau champ : chaque valeur + 5
ecart = 100.0 - champ  # le scalaire se place d'un côté ou de l'autre
carre = champ**2.0  # puissance, valeur par valeur
print(plus_chaud.max("T"), ecart.min("T"), carre.max("q"))  # 35.0 70.0 25.0
print((-champ).min("T"), (2.0**champ).max("q"))  # -30.0 32.0

garde = champ
champ += 1.0  # un nouveau champ, rangé sous le même nom
print(champ.max("T"), garde.max("T"))  # 31.0 30.0 : `garde` a l'ancien

champ += 1.0 construit lui aussi un nouveau champ, et le range sous le même nom : une autre variable qui désignait l’ancien le garde tel quel. Pour modifier un champ lui-même, on passe par les opérations par composante.

Aucune valeur n’est protégée : une division par zéro donne inf, la puissance fractionnaire d’un nombre négatif donne nan, comme avec numpy. La forme pow(champ, s, m), avec un modulo, est refusée.

Sur une composante, en place

add_to_component(c, s), sub_to_component, mul_to_component et div_to_component modifient le champ lui-même : seule la composante c change, sur toutes les zones qui la portent. Elles lèvent une erreur si aucune zone ne porte c, et div_to_component refuse une division par zéro.

champ.add_to_component("T", -1.0)  # en place : T seulement, sur ses deux zones
champ.mul_to_component("q", 2.0)
print(champ.max("T"), champ.max("q"))  # 30.0 12.0

Sur une zone, set_uniform(c, v) donne la valeur v à la composante c en chaque nœud ou point de Gauss.

Entre deux champs

a + b, a - b, a * b, a / b et a ** b combinent les valeurs de deux champs de la même famille. Les valeurs s’associent par support et par composante :

• une composante que les deux champs portent sur un même support se combine valeur à valeur, quelle que soit la zone de ce support qui la range de chaque côté ;
• une composante qu’un seul des deux porte sur un support passe telle quelle, quelle que soit l’opération : a - b rend b, et non -b, sur une composante que seul b porte. C’est ce qui permet de ne toucher qu’à une composante : u + du, où du ne porte que u_x, laisse u_y tel quel, et u * f, où f ne porte que u_x, ne multiplie que u_x ;
• le résultat garde les zones de a, dans l’ordre. Une composante que seul b porte rejoint la zone de a sur ce support quand a n’en a qu’une ; sinon elle garde une zone à elle, comme une zone de b dont le support n’existe pas dans a.

Les deux champs n’ont donc besoin ni des mêmes composantes, ni des mêmes zones :

# a porte T et q, b porte T et p, sur la même zone.
a = pyrucast.NodeField(gauche, ["T", "q"])
a[0].add_to_component("T", 10.0)
a[0].add_to_component("q", 1.0)
b = pyrucast.NodeField(gauche, ["T", "p"])
b[0].add_to_component("T", 3.0)
b[0].add_to_component("p", 4.0)

d = a - b  # T se combine ; q et p, d'un seul côté, passent tels quels
print(d.components())  # ['T', 'q', 'p']
print(d.value(n[0], "T"), d.value(n[0], "q"), d.value(n[0], "p"))  # 7.0 1.0 4.0

champ + zone combine zone, selon les mêmes règles, avec les zones du champ posées sur son support, et recopie les autres ; il lève une erreur si aucune ne l’est. Entre deux zones, les supports doivent être les mêmes, sinon c’est une erreur.

Le résultat respecte les règles de cohérence. Une composante peut toutefois ne pas avoir été obtenue de la même façon d’une zone à l’autre : calculée ici, passée telle quelle là. Deux zones de ce genre qui partagent un nœud n’y donneraient pas la même valeur. Dans ce seul cas, le résultat est vérifié, et l’opération lève une erreur si deux zones ne s’accordent pas en un nœud.

Même support veut dire même objet

Deux zones s’associent quand leurs supports sont le même objet, jamais en comparant leurs nœuds. Deux régions qui se touchent, construites chacune sur sa zone de maillage, ont des supports distincts, même là où elles partagent des nœuds : g + dr laisserait chaque zone telle quelle, et un nœud de la frontière vaudrait 1 dans l’une et 2 dans l’autre. L’opération lève donc une erreur, qui renvoie à restrict_like.

Pour additionner deux régions, on les ramène d’abord sur les mêmes supports. restrict(champ, maillage) pose un champ aux nœuds sur les zones de maillage, et deux champs restreints au même maillage ont les mêmes supports ; restrict_like(b, a) pose b sur les supports de a (voir Opérateurs sur les champs). Un nœud que le champ restreint ne portait pas y vaut 0.

# Deux régions qui se touchent par la diagonale, chacune sur sa zone.
g = pyrucast.NodeField(gauche, ["T"])
g[0].add_to_component("T", 1.0)
dr = pyrucast.NodeField(droite, ["T"])
dr[0].add_to_component("T", 2.0)

try:
    g + dr  # deux supports distincts : n1 vaudrait 1 dans une zone, 2 dans l'autre
except RuntimeError as erreur:
    message = str(erreur)
    print(message)  # incoherent NodeField: node 1, component T: 1 ≠ 2: …

carre = gauche | droite
somme = g.restrict(carre) + dr.restrict(carre)  # les mêmes supports des deux côtés
print(somme.value(n[0], "T"), somme.value(n[1], "T"), somme.value(n[3], "T"))
# 1.0 3.0 2.0 : la diagonale reçoit les deux

Pour un champ aux points de Gauss, le support est le sous-espace : deux espaces construits sur le même maillage sont deux supports (voir Champ aux points de Gauss).

Fonctions mathématiques et masques

pyrucast.field réunit les fonctions qui s’appliquent valeur par valeur et rendent un nouveau champ de la même famille : abs, sqrt, exp, log (népérien), log10, cos, sin, tan, cosh, sinh et tanh. Chacune existe aussi en méthode : champ.sqrt() donne pyrucast.field.sqrt(champ). Elles ne protègent rien non plus : log d’un nombre négatif donne nan.

champ >= s, champ > s, champ <= s et champ < s rendent un masque : un champ de mêmes zones et de mêmes composantes, qui vaut 1 où la condition tient et 0 ailleurs (c’est mask, voir Opérateurs sur les champs). == et != ne comparent pas les valeurs : champ == 0.0 vaut False.

import pyrucast as pc

amorti = pc.field.exp(champ * -0.1)  # e^(-0.1 v), valeur par valeur
print(round(amorti.min("T"), 4))  # 0.0498 : e^-3, au nœud le plus chaud
print(pc.field.cos(champ).max("T"))  # 1.0 : cos 0, au nœud à 0
print(champ.sqrt().max("T"))  # 5.477… : chaque fonction est aussi une méthode

chaud = champ >= 15.0  # 1 où la valeur atteint 15, 0 ailleurs
print(chaud.value(n[0], "T"), chaud.value(n[3], "T"))  # 0.0 1.0

L’union |

L’union compose des zones sans toucher aux valeurs ; c’est l’arithmétique qui combine des valeurs. Elle respecte les règles de cohérence, à sa manière pour chaque famille :

• champ aux nœuds : les zones posées sur un même support sont fondues en une. Un nœud doit avoir la même valeur dans toutes les zones qui le portent, pour chaque composante qu’elles ont en commun, sinon l’union lève une erreur ;
• champ aux points de Gauss : deux zones posées sur un même sous-espace restent côte à côte si leurs composantes sont différentes ; une même composante portée deux fois sur un sous-espace est refusée.

consolidate() fond, à la demande, les zones qui partagent un support (voir Champ aux nœuds et Champ aux points de Gauss).

Ce qui est mis en cache

Un champ ne garde aucun cache : une zone ne contient que le nom de ses composantes, ses valeurs et une référence à son support. Ce qui ne se calcule qu’une fois vit sur les supports, et c’est ce qui rend les opérations ci-dessus possibles ou rapides :

• le nuage de nœuds d’une zone de maillage. NodeField(maillage, …), restrict et les opérateurs qui produisent un champ aux nœuds posent leurs zones sur le nuage des nœuds de chaque zone de maillage. Ce nuage est calculé une fois, puis rendu tel quel à chaque demande, tant que la zone de maillage ne change pas. Deux champs construits sur le même maillage reçoivent donc le même support et se combinent nœud à nœud : c’est ce qui fait marcher a - b et restrict(…, carre) dans les exemples ci-dessus ;
• la table qui donne la ligne d’un nœud dans ce nuage, construite à la première lecture d’un nœud : value(nœud, c), set_value et les opérateurs trouvent la ligne d’un nœud sans parcourir la liste ;
• les tables de l’espace éléments finis (points et poids de Gauss, fonctions de forme et leurs dérivées) et la coloration de ses mailles, pour un champ aux points de Gauss (voir Champ aux points de Gauss).

Côté Rust

Le socle commun est écrit une seule fois, dans deux traits de src/containers/field.rs. Python n’en voit que le résultat : les méthodes et les opérateurs décrits plus haut, présents sur les quatre classes.

• SubField, la zone. Un type n’a que cinq choses à fournir : support(), components(), values(), values_mut() et same_support_with(…), une zone neuve, à zéro, sur le même support. Tout le reste en découle : component_index, min, max, sum, xtx, dot, pscal, map_all, add_to_component…, merge_components, select_components, rename_component. SubNodeField et SubElementField l’implémentent. SubNodeField fournit en plus row_points(), le nœud de chaque ligne, lu sur place dans le support : c’est ce qui permet de compter un nœud une fois (la première zone qui le porte le compte, les suivantes le sautent) et de vérifier la deuxième règle. Une zone aux points de Gauss garde la réponse par défaut, None.
• Field, l’agrégat. Il est implémenté d’office pour tout Aggregate dont la zone est un SubField : NodeField et ElementField le reçoivent sans une ligne de code : components, min, max, sum, map_subs, combine_scalar, map_all, merge_field, merge_subfield, dot_field, pscal_field, filter_components, rename_component. Les réductions comptent une fois un nœud que plusieurs zones portent ; merge_field, merge_subfield, dot_field et pscal_field associent les valeurs par support et par composante.
• Règles de cohérence. add_sub refuse une zone qui porterait une composante déjà portée sur son support : c’est le crochet check_add, qu’une union ne passe pas, puisque son finalize vérifie ou fond le résultat entier. Les opérateurs qui empilent des contributions partielles avant de les additionner nœud à nœud (divergence, forces internes et externes, sum_to_nodes) ne vérifient pas la deuxième règle en chemin : ce ne sont pas encore des champs.

Le contrat ne dit pas ce qu’est une ligne — un nœud pour SubNodeField, un couple (maille, point de Gauss) pour SubElementField. Il fixe seulement la disposition des valeurs dans un tableau plat, où la composante varie le plus vite :

   zone aux composantes [c0, c1]
   values = [ ligne0.c0, ligne0.c1, ligne1.c0, ligne1.c1, … ]

#[test]
fn le_contrat_des_champs_s_ecrit_une_fois() -> Result<()> {
    let coords = Handle::new(Coords::new(1)?);
    let a = Node::create_in(coords.clone(), &[0.0])?;
    let b = Node::create_in(coords, &[1.0])?;
    let support = mesh::poi1_from_nodes(&[a.clone(), b.clone()])?;

    let t = NodeField::new(&support, vec!["T".into()])?;
    t.add_to_component("T", 2.0)?; // Field : en place, sur chaque zone qui porte T
    let tq = NodeField::new(&support, vec!["T".into(), "q".into()])?;
    tq.add_to_component("T", 3.0)?;

    // Les opérateurs rendent un `Result` : combiner deux champs peut échouer.
    let s = (&t + &tq)?; // même support : T se combine, q passe tel quel
    assert_eq!(s.value(a.id(), "T")?, 5.0);
    assert_eq!((2.0 * &t)?.value(b.id(), "T")?, 4.0);

    // Rust n'a pas d'opérateur de puissance : on passe par les primitives.
    let carre = t.combine_scalar(f64::powf, 2.0)?;
    let puissance = t.merge_field(&tq, f64::powf)?; // 2³
    assert_eq!(carre.value(a.id(), "T")?, 4.0);
    assert_eq!(puissance.value(a.id(), "T")?, 8.0);

    // Une zone se lit par le trait `SubField`.
    let zone = t.get(0)?;
    assert_eq!(zone.read().component_index("T"), Some(0));
    assert_eq!(SubField::max(&*zone.read(), Some("T"))?, 2.0);
    Ok(())
}

• Opérateurs. Trois macros écrivent + - * / : impl_subfield_scalar_ops! (une zone et un scalaire, d’un côté ou de l’autre), impl_subfield_field_ops! (deux zones, par merge_components) et impl_field_ops! (deux agrégats par merge_field, un agrégat et un scalaire par combine_scalar). Combiner peut échouer, sauf une zone avec un scalaire : ces opérateurs rendent un Result. Une zone prise par valeur est modifiée sur place, sans copie ; prise par référence, elle est d’abord copiée.
• Puissance. Rust n’a pas d’opérateur ** : combine_scalar(f64::powf, s) et merge_field(&autre, f64::powf) font ce que fait ** en Python.
• Support. same_support compare les Handle par identité (Handle::same_object).
• Mêmes composantes exigées. check_same_components refuse deux zones dont les supports ou les composantes diffèrent. L’interpolation d’une Evolution l’appelle avant merge_components : deux valeurs tabulées d’un même champ doivent porter les mêmes composantes.
• Lecture par rang. Un noyau de calcul ne cherche pas une composante par son nom à chaque point : resolve_components traduit, une fois par zone, une liste de noms en rangs dans la zone, et resolve_optional_components marque ABSENT_COMPONENT les noms que la zone ne porte pas.
• Fonctions mathématiques. Les traits MapValues et Pscal donnent un nom commun aux quatre types (zone ou agrégat, nœuds ou points de Gauss) : ops::field::cos, ops::field::psca… sont écrits une seule fois.
• Vue de lecture. Field::view prend un verrou de lecture par zone et le garde le temps de la vue. Les opérateurs qui lisent beaucoup (gradient, deformation, divergence, restrict, integral, le tracé…) lisent ainsi les valeurs sur place, sans copie ; d’autres lectures peuvent avoir lieu en même temps, une écriture attend la fin de la vue.
• Parallélisme. min, max, sum, xtx, map_all… parcourent les valeurs en parallèle. min et max donnent le même résultat quel que soit le nombre de threads ; sum, xtx et dot peuvent différer au dernier chiffre.