64#include "utilitaires.h"
67#include "param_elliptic.h"
70#include "sym_tensor.h"
94 if (
aaa.get_etat() == ETATZERO) {
109 int nz =
mgrid.get_nzone() ;
124 int nt =
mgrid.get_nt(0) ;
125 int np =
mgrid.get_np(0) ;
132 source.set_spectral_va().ylm() ;
139 tilde_mu.set_spectral_va().set_etat_cf_qcq() ;
140 tilde_mu.set_spectral_va().c_cf->annule_hard() ;
141 x_new.annule_hard() ;
142 x_new.set_spectral_base(base) ;
143 x_new.set_spectral_va().set_etat_cf_qcq() ;
144 x_new.set_spectral_va().c_cf->annule_hard() ;
167 for (
int k=0 ;
k<np+1 ;
k++) {
168 for (
int j=0 ;
j<nt ;
j++) {
171 if ( (nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
176 for (
int i=0;
i<nr;
i++) {
180 for (
int i=0;
i<nr;
i++) {
201 for (
int k=0 ;
k<np+1 ;
k++) {
202 for (
int j=0 ;
j<nt ;
j++) {
205 if ( (nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
245 for (
int i=0;
i<nr;
i++) {
252 for (
int i=0;
i<nr;
i++) {
259 for (
int i=0;
i<nr;
i++) {
279 for (
int k=0 ;
k<np+1 ;
k++) {
280 for (
int j=0 ;
j<nt ;
j++) {
283 if ( (nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
327 for (
int i=0;
i<nr;
i++) {
334 for (
int i=0;
i<nr;
i++) {
346 int taille = 2*
nz_bc ;
347 if (
cedbc) taille-- ;
378 for (
int k=0 ;
k<np+1 ;
k++)
379 for (
int j=0 ;
j<nt ;
j++) {
381 if ((nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
389 int nr =
mgrid.get_nr(1) ;
420 double blob = (*boundmu).val_in_bound_jk(1,
j,
k);
542 nr =
mgrid.get_nr(0) ;
543 for (
int i=0 ;
i<nr ;
i++) {
545 mw.set(0,
k,
j,
i) = 0 ;
547 nr =
mgrid.get_nr(1) ;
548 for (
int i=0 ;
i<nr ;
i++) {
559 for (
int i=0 ;
i<nr ;
i++) {
572 for (
int i=0 ;
i<nr ;
i++) {
583 if (
tilde_mu.set_spectral_va().c != 0x0)
584 delete tilde_mu.set_spectral_va().c ;
585 tilde_mu.set_spectral_va().c = 0x0 ;
586 tilde_mu.set_spectral_va().ylm_i() ;
588 if (
x_new.set_spectral_va().c != 0x0)
589 delete x_new.set_spectral_va().c ;
590 x_new.set_spectral_va().c = 0x0 ;
591 x_new.set_spectral_va().ylm_i();
612 int nz =
mgrid.get_nzone() ;
618 if (
par_bc->get_n_int() > 0)
628 int nt =
mgrid.get_nt(0) ;
629 int np =
mgrid.get_np(0) ;
670 if ( (
bbt.get_etat() == ETATZERO) && (hh.
get_etat() == ETATZERO) ) {
689 source.set_spectral_va().ylm() ;
691 bool bnull = (
bbt.get_etat() == ETATZERO) ;
696 hoverr.set_spectral_va().ylm() ;
698 dhdr.set_spectral_va().ylm() ;
700 h_coq.set_spectral_va().ylm() ;
703 dh_coq.set_spectral_va().ylm() ;
714 if ((
par_mat->get_n_int_mod() >= 4)
715 &&(
par_mat->get_n_tbl_mod()>=1)
716 &&(
par_mat->get_n_matrice_mod()>=1)
717 &&(
par_mat->get_n_itbl_mod()>=1)) {
750 pnr->set_etat_qcq() ;
752 for (
int l=0;
l<nz;
l++)
766 hrr.set_spectral_base(base) ;
767 hrr.set_spectral_va().set_etat_cf_qcq() ;
768 hrr.set_spectral_va().c_cf->annule_hard() ;
771 tilde_eta.set_spectral_va().set_etat_cf_qcq() ;
772 tilde_eta.set_spectral_va().c_cf->annule_hard() ;
810 for (
int k=0 ;
k<np+1 ;
k++) {
811 for (
int j=0 ;
j<nt ;
j++) {
814 if ( (nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
855 ope.set(2*nr-1,
col) = 0 ;
856 ope.set(3*nr-1,
col) = 0 ;
868 ope.set(nr-1, nr-1) = 1 ;
869 ope.set(2*nr-1, 2*nr-1) = 1 ;
870 ope.set(3*nr-1, 3*nr-1) = 1 ;
915 sec.set(3*nr-1) = 0 ;
917 for (
int i=0;
i<nr;
i++) {
933 for (
int i=0;
i<nr;
i++) {
958 for (
int k=0 ;
k<np+1 ;
k++) {
959 for (
int j=0 ;
j<nt ;
j++) {
962 if ( (nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
1021 sec.set_etat_qcq() ;
1033 sec.set(
lin+nr) = -0.5*(*
h_coq.get_spectral_va().c_cf)
1053 for (
int i=0;
i<nr;
i++) {
1061 for (
int i=0;
i<nr;
i++) {
1069 for (
int i=0;
i<nr;
i++) {
1077 for (
int i=0;
i<nr;
i++) {
1099 for (
int k=0 ;
k<np+1 ;
k++) {
1100 for (
int j=0 ;
j<nt ;
j++) {
1103 if ( (nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
1105 ope.set_etat_qcq() ;
1165 sec.set_etat_qcq() ;
1177 sec.set(
lin+nr) = -0.5*(*
hoverr.get_spectral_va().c_cf)
1199 for (
int i=0;
i<nr;
i++) {
1207 for (
int i=0;
i<nr;
i++) {
1215 for (
int i=0;
i<nr;
i++) {
1229 int taille = 3*
nz_bc ;
1230 if (
cedbc) taille-- ;
1276 for (
int k=0 ;
k<np+1 ;
k++)
1277 for (
int j=0 ;
j<nt ;
j++) {
1279 if ((nullite_plm(
j, nt,
k, np, base) == 1) && (
l_q > 1)) {
1287 int nr =
mgrid.get_nr(1);
1295 double blob = (*(
bound_hrr.get_spectral_va().c_cf)).val_in_bound_jk(1,
j,
k);
1296 double blob2 = (*(
bound_eta.get_spectral_va().c_cf)).val_in_bound_jk(1,
j,
k);
1513 nr =
mgrid.get_nr(0) ;
1514 for (
int i=0 ;
i<nr ;
i++) {
1517 mw.set(0,
k,
j,
i) = 0 ;
1521 for (
int i=0 ;
i<nr ;
i++) {
1541 for (
int i=0 ;
i<nr ;
i++) {
1559 if (
hrr.set_spectral_va().c != 0x0)
1560 delete hrr.set_spectral_va().c;
1561 hrr.set_spectral_va().c = 0x0 ;
1562 hrr.set_spectral_va().ylm_i() ;
1564 if (
tilde_eta.set_spectral_va().c != 0x0)
1595 int nz =
mgrid.get_nzone() ;
1605 int nt = mgrid.
get_nt(0) ;
1606 int np = mgrid.
get_np(0) ;
1608 Scalar source = hh ;
1610 Scalar source_coq = hh ;
1612 source_coq.set_spectral_va().ylm() ;
1613 Base_val base = source.get_spectral_base() ;
1615 int lmax = base.give_lmax(mgrid, 0) + 1;
1622 Mtbl_cf sol_part_hrr(mgrid, base) ; sol_part_hrr.annule_hard() ;
1623 Mtbl_cf sol_part_eta(mgrid, base) ; sol_part_eta.annule_hard() ;
1624 Mtbl_cf sol_hom1_hrr(mgrid, base) ; sol_hom1_hrr.annule_hard() ;
1625 Mtbl_cf sol_hom1_eta(mgrid, base) ; sol_hom1_eta.annule_hard() ;
1626 Mtbl_cf sol_hom2_hrr(mgrid, base) ; sol_hom2_hrr.annule_hard() ;
1627 Mtbl_cf sol_hom2_eta(mgrid, base) ; sol_hom2_eta.annule_hard() ;
1629 bool need_calculation = true ;
1634 int l_q, m_q, base_r ;
1635 Itbl mat_done(lmax) ;
1641 int nr = mgrid.
get_nr(lz) ;
1642 double alpha = mp_aff->
get_alpha()[lz] ;
1643 Matrice ope(2*nr, 2*nr) ;
1644 if (need_calculation && (par_mat != 0x0)) mat_done.annule_hard() ;
1646 for (
int k=0 ; k<np+1 ; k++) {
1647 for (
int j=0 ; j<nt ; j++) {
1649 base.give_quant_numbers(lz, k, j, m_q, l_q, base_r) ;
1650 if ( (nullite_plm(j, nt, k, np, base) == 1) && (l_q < 2)) {
1651 if (need_calculation) {
1652 ope.set_etat_qcq() ;
1653 Diff_dsdx od(base_r, nr) ;
const Matrice& md = od.get_matrice() ;
1654 Diff_sx os(base_r, nr) ;
const Matrice& ms = os.get_matrice() ;
1656 for (
int lin=0; lin<nr; lin++)
1657 for (
int col=0; col<nr; col++)
1658 ope.set(lin,col) = md(lin,col) + 3*ms(lin,col) ;
1659 for (
int lin=0; lin<nr; lin++)
1660 for (
int col=0; col<nr; col++)
1661 ope.set(lin,col+nr) = -l_q*(l_q+1)*ms(lin,col) ;
1662 for (
int lin=0; lin<nr; lin++)
1663 for (
int col=0; col<nr; col++)
1664 ope.set(lin+nr,col) = -0.5*ms(lin,col) ;
1665 for (
int lin=0; lin<nr; lin++)
1666 for (
int col=0; col<nr; col++)
1667 ope.set(lin+nr,col+nr) = md(lin,col) + 3*ms(lin, col);
1670 for (
int col=0; col<2*nr; col++) {
1671 ope.set(nr-1, col) = 0 ;
1672 ope.set(2*nr-1, col) = 0 ;
1676 for (
int col=0; col<nr; col++) {
1677 ope.set(nr-1, col) = pari ;
1678 ope.set(2*nr-1, col+nr) = pari ;
1683 ope.set(nr-1, nr-1) = 1 ;
1684 ope.set(2*nr-1, 2*nr-1) = 1 ;
1688 if ((par_mat != 0x0) && (mat_done(l_q) == 0)) {
1689 Matrice* pope =
new Matrice(ope) ;
1691 mat_done.set(l_q) = 1 ;
1695 const Matrice& oper = (par_mat == 0x0 ? ope :
1698 sec.set_etat_qcq() ;
1699 for (
int lin=0; lin<nr; lin++)
1700 sec.set(lin) = (*source.get_spectral_va().c_cf)(lz, k, j, lin) ;
1701 for (
int lin=0; lin<nr; lin++)
1702 sec.set(nr+lin) = -0.5*(*source.get_spectral_va().c_cf)
1708 for (
int col=0; col<nr; col++) {
1710 (*source_coq.get_spectral_va().c_cf)(lz, k, j, col) ;
1713 sec.set(nr-1) = h0 / 3. ;
1715 sec.set(2*nr-1) = 0 ;
1716 Tbl sol = oper.inverse(sec) ;
1717 for (
int i=0; i<nr; i++) {
1718 sol_part_hrr.set(lz, k, j, i) = sol(i) ;
1719 sol_part_eta.set(lz, k, j, i) = sol(i+nr) ;
1731 for (
int lz=1; lz<nz-1; lz++) {
1732 if (need_calculation && (par_mat != 0x0)) mat_done.annule_hard() ;
1733 int nr = mgrid.
get_nr(lz) ;
1736 assert(mgrid.
get_nt(lz) == nt) ;
1737 assert(mgrid.
get_np(lz) == np) ;
1738 double alpha = mp_aff->
get_alpha()[lz] ;
1739 double ech = mp_aff->
get_beta()[lz] / alpha ;
1740 Matrice ope(2*nr, 2*nr) ;
1742 for (
int k=0 ; k<np+1 ; k++) {
1743 for (
int j=0 ; j<nt ; j++) {
1745 base.give_quant_numbers(lz, k, j, m_q, l_q, base_r) ;
1746 if ( (nullite_plm(j, nt, k, np, base) == 1) && (l_q < 2)) {
1747 if (need_calculation) {
1748 ope.set_etat_qcq() ;
1749 Diff_xdsdx oxd(base_r, nr) ;
const Matrice& mxd = oxd.get_matrice() ;
1750 Diff_dsdx od(base_r, nr) ;
const Matrice& md = od.get_matrice() ;
1751 Diff_id oid(base_r, nr) ;
const Matrice& mid = oid.get_matrice() ;
1753 for (
int lin=0; lin<nr; lin++)
1754 for (
int col=0; col<nr; col++)
1755 ope.set(lin,col) = mxd(lin,col) + ech*md(lin,col)
1757 for (
int lin=0; lin<nr; lin++)
1758 for (
int col=0; col<nr; col++)
1759 ope.set(lin,col+nr) = -l_q*(l_q+1)*mid(lin,col) ;
1760 for (
int lin=0; lin<nr; lin++)
1761 for (
int col=0; col<nr; col++)
1762 ope.set(lin+nr,col) = -0.5*mid(lin,col) ;
1763 for (
int lin=0; lin<nr; lin++)
1764 for (
int col=0; col<nr; col++)
1765 ope.set(lin+nr,col+nr) = mxd(lin,col) + ech*md(lin,col)
1768 for (
int col=0; col<2*nr; col++) {
1769 ope.set(ind0+nr-1, col) = 0 ;
1770 ope.set(ind1+nr-1, col) = 0 ;
1772 ope.set(ind0+nr-1, ind0) = 1 ;
1773 ope.set(ind1+nr-1, ind1) = 1 ;
1776 if ((par_mat != 0x0) && (mat_done(l_q) == 0)) {
1777 Matrice* pope =
new Matrice(ope) ;
1779 mat_done.set(l_q) = 1 ;
1782 const Matrice& oper = (par_mat == 0x0 ? ope :
1785 sec.set_etat_qcq() ;
1786 for (
int lin=0; lin<nr; lin++)
1787 sec.set(lin) = (*source_coq.get_spectral_va().c_cf)
1789 for (
int lin=0; lin<nr; lin++)
1790 sec.set(nr+lin) = -0.5*(*source_coq.get_spectral_va().c_cf)
1792 sec.set(ind0+nr-1) = 0 ;
1793 sec.set(ind1+nr-1) = 0 ;
1794 Tbl sol = oper.inverse(sec) ;
1796 for (
int i=0; i<nr; i++) {
1797 sol_part_hrr.set(lz, k, j, i) = sol(i) ;
1798 sol_part_eta.set(lz, k, j, i) = sol(i+nr) ;
1801 sec.set(ind0+nr-1) = 1 ;
1802 sol = oper.inverse(sec) ;
1803 for (
int i=0; i<nr; i++) {
1804 sol_hom1_hrr.set(lz, k, j, i) = sol(i) ;
1805 sol_hom1_eta.set(lz, k, j, i) = sol(i+nr) ;
1807 sec.set(ind0+nr-1) = 0 ;
1808 sec.set(ind1+nr-1) = 1 ;
1809 sol = oper.inverse(sec) ;
1810 for (
int i=0; i<nr; i++) {
1811 sol_hom2_hrr.set(lz, k, j, i) = sol(i) ;
1812 sol_hom2_eta.set(lz, k, j, i) = sol(i+nr) ;
1823 if (need_calculation && (par_mat != 0x0)) mat_done.annule_hard() ;
1824 int nr = mgrid.
get_nr(lz) ;
1827 assert(mgrid.
get_nt(lz) == nt) ;
1828 assert(mgrid.
get_np(lz) == np) ;
1829 double alpha = mp_aff->
get_alpha()[lz] ;
1830 Matrice ope(2*nr, 2*nr) ;
1832 for (
int k=0 ; k<np+1 ; k++) {
1833 for (
int j=0 ; j<nt ; j++) {
1835 base.give_quant_numbers(lz, k, j, m_q, l_q, base_r) ;
1836 if ( (nullite_plm(j, nt, k, np, base) == 1) && (l_q < 2)) {
1837 if (need_calculation) {
1838 ope.set_etat_qcq() ;
1839 Diff_dsdx od(base_r, nr) ;
const Matrice& md = od.get_matrice() ;
1840 Diff_sx os(base_r, nr) ;
const Matrice& ms = os.get_matrice() ;
1842 for (
int lin=0; lin<nr; lin++)
1843 for (
int col=0; col<nr; col++)
1844 ope.set(lin,col) = - md(lin,col) + 3*ms(lin,col) ;
1845 for (
int lin=0; lin<nr; lin++)
1846 for (
int col=0; col<nr; col++)
1847 ope.set(lin,col+nr) = -l_q*(l_q+1)*ms(lin,col) ;
1848 for (
int lin=0; lin<nr; lin++)
1849 for (
int col=0; col<nr; col++)
1850 ope.set(lin+nr,col) = -0.5*ms(lin,col) ;
1851 for (
int lin=0; lin<nr; lin++)
1852 for (
int col=0; col<nr; col++)
1853 ope.set(lin+nr,col+nr) = -md(lin,col) + 3*ms(lin, col) ;
1856 for (
int col=0; col<2*nr; col++) {
1857 ope.set(ind0+nr-2, col) = 0 ;
1858 ope.set(ind0+nr-1, col) = 0 ;
1859 ope.set(ind1+nr-2, col) = 0 ;
1860 ope.set(ind1+nr-1, col) = 0 ;
1862 for (
int col=0; col<nr; col++) {
1863 ope.set(ind0+nr-1, col+ind0) = 1 ;
1864 ope.set(ind1+nr-1, col+ind1) = 1 ;
1866 ope.set(ind0+nr-2, ind0+1) = 1 ;
1867 ope.set(ind1+nr-2, ind1+1) = 1 ;
1870 if ((par_mat != 0x0) && (mat_done(l_q) == 0)) {
1871 Matrice* pope =
new Matrice(ope) ;
1873 mat_done.set(l_q) = 1 ;
1876 const Matrice& oper = (par_mat == 0x0 ? ope :
1879 sec.set_etat_qcq() ;
1880 for (
int lin=0; lin<nr; lin++)
1881 sec.set(lin) = (*source.get_spectral_va().c_cf)
1883 for (
int lin=0; lin<nr; lin++)
1884 sec.set(nr+lin) = -0.5*(*source.get_spectral_va().c_cf)
1886 sec.set(ind0+nr-2) = 0 ;
1887 sec.set(ind0+nr-1) = 0 ;
1888 sec.set(ind1+nr-2) = 0 ;
1889 sec.set(ind1+nr-1) = 0 ;
1890 Tbl sol = oper.inverse(sec) ;
1891 for (
int i=0; i<nr; i++) {
1892 sol_part_hrr.set(lz, k, j, i) = sol(i) ;
1893 sol_part_eta.set(lz, k, j, i) = sol(i+nr) ;
1896 sec.set(ind0+nr-2) = 1 ;
1897 sol = oper.inverse(sec) ;
1898 for (
int i=0; i<nr; i++) {
1899 sol_hom1_hrr.set(lz, k, j, i) = sol(i) ;
1900 sol_hom1_eta.set(lz, k, j, i) = sol(i+nr) ;
1902 sec.set(ind0+nr-2) = 0 ;
1903 sec.set(ind1+nr-2) = 1 ;
1904 sol = oper.inverse(sec) ;
1905 for (
int i=0; i<nr; i++) {
1906 sol_hom2_hrr.set(lz, k, j, i) = sol(i) ;
1907 sol_hom2_eta.set(lz, k, j, i) = sol(i+nr) ;
1916 Mtbl_cf dhom1_hrr = sol_hom1_hrr ;
1917 Mtbl_cf dhom2_hrr = sol_hom2_hrr ;
1918 Mtbl_cf dpart_hrr = sol_part_hrr ;
1919 Mtbl_cf dhom1_eta = sol_hom1_eta ;
1920 Mtbl_cf dhom2_eta = sol_hom2_eta ;
1921 Mtbl_cf dpart_eta = sol_part_eta ;
1923 dhom1_hrr.dsdx() ; dhom1_eta.dsdx() ;
1924 dhom2_hrr.dsdx() ; dhom2_eta.dsdx() ;
1925 dpart_hrr.dsdx() ; dpart_eta.dsdx() ;
1929 int taille = 2*(nz-1) ;
1933 Tbl sec_membre(taille) ;
1934 Matrice systeme(taille, taille) ;
1935 int ligne ;
int colonne ;
1939 for (
int k=0 ; k<np+1 ; k++)
1940 for (
int j=0 ; j<nt ; j++) {
1941 base.give_quant_numbers(0, k, j, m_q, l_q, base_r) ;
1942 if ((nullite_plm(j, nt, k, np, base) == 1) && (l_q < 2)) {
1945 systeme.annule_hard() ;
1946 sec_membre.annule_hard() ;
1953 sec_membre.set(ligne) = 0.;
1956 sec_membre.set(ligne) = 0.;
1961 double alpha2 = mp_aff->
get_alpha()[1] ;
1962 double blob = (*(bound_hrr.
get_spectral_va().c_cf)).val_in_bound_jk(1, j, k);
1963 double blob2 = (*(bound_eta.
get_spectral_va().c_cf)).val_in_bound_jk(1, j, k);
1968 systeme.set(ligne, colonne) = - (4. - double(l_q*(l_q+1.)))*sol_hom1_hrr.val_in_bound_jk(1,j,k) - 2.*dhom1_hrr.val_in_bound_jk(1,j,k)/alpha2;
1970 systeme.set(ligne, colonne +1) = - (4. - double(l_q*(l_q+1.)))*sol_hom2_hrr.val_in_bound_jk(1,j,k) - 2.*dhom2_hrr.val_in_bound_jk(1,j,k)/alpha2;
1972 sec_membre.set(ligne) = (4. - double(l_q*(l_q+1.)))*sol_part_hrr.val_in_bound_jk(1,j,k) + 2.*dpart_hrr.val_in_bound_jk(1,j,k)/alpha2 - blob;
1978 systeme.set(ligne, colonne) = - ( 2.*dhom1_eta.val_in_bound_jk(1,j,k)/alpha2 + (2.-double(l_q*(l_q+1.)))*sol_hom1_eta.val_in_bound_jk(1,j,k) + 2.*sol_hom1_hrr.val_in_bound_jk(1,j,k));
1980 systeme.set(ligne, colonne +1) = - ( 2.*dhom2_eta.val_in_bound_jk(1,j,k)/alpha2 + (2.-double(l_q*(l_q+1.)))*sol_hom2_eta.val_in_bound_jk(1,j,k) + 2.*sol_hom2_hrr.val_in_bound_jk(1,j,k));
1982 sec_membre.set(ligne) = 2.*dpart_eta.val_in_bound_jk(1,j,k)/alpha2 + (2.-double(l_q*(l_q+1.)))*sol_part_eta.val_in_bound_jk(1,j,k) + 2.*sol_part_hrr.val_in_bound_jk(1,j,k) - blob2;
1987 systeme.set(ligne, colonne) =
1988 sol_hom1_hrr.val_out_bound_jk(1, j, k) ;
1989 systeme.set(ligne, colonne+1) =
1990 sol_hom2_hrr.val_out_bound_jk(1, j, k) ;
1992 sec_membre.set(ligne) -= sol_part_hrr.val_out_bound_jk(1, j, k) ;
1995 systeme.set(ligne, colonne) =
1996 sol_hom1_eta.val_out_bound_jk(1, j, k) ;
1997 systeme.set(ligne, colonne+1) =
1998 sol_hom2_eta.val_out_bound_jk(1, j, k) ;
2000 sec_membre.set(ligne) -= sol_part_eta.val_out_bound_jk(1, j, k) ;
2006 for (
int zone=2 ; zone<nz-1 ; zone++) {
2007 nr = mgrid.
get_nr(zone) ;
2011 systeme.set(ligne, colonne) =
2012 - sol_hom1_hrr.val_in_bound_jk(zone, j, k) ;
2013 systeme.set(ligne, colonne+1) =
2014 - sol_hom2_hrr.val_in_bound_jk(zone, j, k) ;
2016 sec_membre.set(ligne) += sol_part_hrr.val_in_bound_jk(zone, j, k) ;
2019 systeme.set(ligne, colonne) =
2020 - sol_hom1_eta.val_in_bound_jk(zone, j, k) ;
2021 systeme.set(ligne, colonne+1) =
2022 - sol_hom2_eta.val_in_bound_jk(zone, j, k) ;
2024 sec_membre.set(ligne) += sol_part_eta.val_in_bound_jk(zone, j, k) ;
2028 systeme.set(ligne, colonne) =
2029 sol_hom1_hrr.val_out_bound_jk(zone, j, k) ;
2030 systeme.set(ligne, colonne+1) =
2031 sol_hom2_hrr.val_out_bound_jk(zone, j, k) ;
2033 sec_membre.set(ligne) -= sol_part_hrr.val_out_bound_jk(zone, j, k) ;
2036 systeme.set(ligne, colonne) =
2037 sol_hom1_eta.val_out_bound_jk(zone, j, k) ;
2038 systeme.set(ligne, colonne+1) =
2039 sol_hom2_eta.val_out_bound_jk(zone, j, k) ;
2041 sec_membre.set(ligne) -= sol_part_eta.val_out_bound_jk(zone, j, k) ;
2047 nr = mgrid.
get_nr(nz-1) ;
2051 systeme.set(ligne, colonne) =
2052 - sol_hom1_hrr.val_in_bound_jk(nz-1, j, k) ;
2053 systeme.set(ligne, colonne+1) =
2054 - sol_hom2_hrr.val_in_bound_jk(nz-1, j, k) ;
2056 sec_membre.set(ligne) += sol_part_hrr.val_in_bound_jk(nz-1, j, k) ;
2059 systeme.set(ligne, colonne) =
2060 - sol_hom1_eta.val_in_bound_jk(nz-1, j, k) ;
2061 systeme.set(ligne, colonne+1) =
2062 - sol_hom2_eta.val_in_bound_jk(nz-1, j, k) ;
2064 sec_membre.set(ligne) += sol_part_eta.val_in_bound_jk(nz-1, j, k) ;
2070 Tbl facteur = systeme.inverse(sec_membre) ;
2076 for (
int i=0 ; i<nr ; i++) {
2077 mhrr.set(0, k, j, i) = sol_part_hrr(0, k, j, i) ;
2078 meta.set(0, k, j, i) = sol_part_eta(0, k, j, i) ;
2080 for (
int zone=1 ; zone<nz-1 ; zone++) {
2081 nr = mgrid.
get_nr(zone) ;
2082 for (
int i=0 ; i<nr ; i++) {
2083 mhrr.set(zone, k, j, i) = sol_part_hrr(zone, k, j, i)
2084 + facteur(conte)*sol_hom1_hrr(zone, k, j, i)
2085 + facteur(conte+1)*sol_hom2_hrr(zone, k, j, i) ;
2087 meta.set(zone, k, j, i) = sol_part_eta(zone, k, j, i)
2088 + facteur(conte)*sol_hom1_eta(zone, k, j, i)
2089 + facteur(conte+1)*sol_hom2_eta(zone, k, j, i) ;
2093 nr = mgrid.
get_nr(nz-1) ;
2094 for (
int i=0 ; i<nr ; i++) {
2095 mhrr.set(nz-1, k, j, i) = sol_part_hrr(nz-1, k, j, i)
2096 + facteur(conte)*sol_hom1_hrr(nz-1, k, j, i)
2097 + facteur(conte+1)*sol_hom2_hrr(nz-1, k, j, i) ;
2099 meta.set(nz-1, k, j, i) = sol_part_eta(nz-1, k, j, i)
2100 + facteur(conte)*sol_hom1_eta(nz-1, k, j, i)
2101 + facteur(conte+1)*sol_hom2_eta(nz-1, k, j, i) ;
Bases of the spectral expansions.
void mult_x()
The basis is transformed as with a multiplication by .
int give_lmax(const Mg3d &mgrid, int lz) const
Returns the highest multipole for a given grid.
void give_quant_numbers(int, int, int, int &, int &, int &) const
Computes the various quantum numbers and 1d radial base.
Class for the elementary differential operator (see the base class Diff ).
Class for the elementary differential operator Identity (see the base class Diff ).
Class for the elementary differential operator division by (see the base class Diff ).
Class for the elementary differential operator (see the base class Diff ).
Time evolution with partial storage (*** under development ***).
Basic integer array class.
const double * get_beta() const
Returns the pointer on the array beta.
const double * get_alpha() const
Returns the pointer on the array alpha.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Coefficients storage for the multi-domain spectral method.
void annule_hard()
Sets the Mtbl_cf to zero in a hard way.
int get_n_matrice_mod() const
Returns the number of modifiable Matrice 's addresses in the list.
void add_matrice_mod(Matrice &ti, int position=0)
Adds the address of a new modifiable Matrice to the list.
Matrice & get_matrice_mod(int position=0) const
Returns the reference of a modifiable Matrice stored in the list.
Tensor field of valence 0 (or component of a tensorial field).
void div_r_dzpuis(int ced_mult_r)
Division by r everywhere but with the output flag dzpuis set to ced_mult_r .
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
const Scalar & dsdr() const
Returns of *this .
Valeur & set_spectral_va()
Returns va (read/write version)
const Valeur & get_spectral_va() const
Returns va (read only version)
void annule_hard()
Sets the Scalar to zero in a hard way.
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
const Base_val & get_spectral_base() const
Returns the spectral bases of the Valeur va
void set_spectral_base(const Base_val &)
Sets the spectral bases of the Valeur va
void sol_Dirac_A2(const Scalar &aaa, Scalar &tilde_mu, Scalar &x_new, Scalar bound_mu, const Param *par_bc)
Same resolution as sol_Dirac_Abound, but here the boundary conditions are the degenerate elliptic con...
void sol_Dirac_BC3(const Scalar &bb, const Scalar &hh, Scalar &hrr, Scalar &tilde_eta, Scalar &ww, Scalar bound_hrr, Scalar bound_eta, Param *par_bc, Param *par_mat)
Same resolution as sol_Dirac_Abound, but here the boundary conditions are the degenerate elliptic con...
void set_etat_cf_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl_cf c_c...
void ylm()
Computes the coefficients of *this.
Mtbl * c
Values of the function at the points of the multi-grid
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
void ylm_i()
Inverse of ylm()
#define R_CHEBP
base de Cheb. paire (rare) seulement
const Map *const mp
Mapping on which the numerical values at the grid points are defined.