docs: def write_STELLA_model(self,name):
          """
          Write an initial model in a format that may easily be read by the
          radiation hydrodynamics code STELLA.
  
          Parameters
          ----------
          name : string
              an identifier for the model. There are two output files from
              this method, which will be <name>.hyd and <name>.abn, which
              contain the profiles for the hydro and abundance variables,
              respectively.
  
          """
  
          # Hydro variables:
          zn = np.array(self.get('zone'),np.int64)
          Mr = self.get('mass')[::-1]
          dM = 10. ** self.get('logdq')[::-1] * self.header_attr['star_mass']
          R = self.get('radius')[::-1] * ast.rsun_cm
          dR = np.insert( np.diff(R), 0, R[0] )
          Rho = 10. ** self.get('logRho')[::-1]
          PRE = 10. ** self.get('logP')[::-1]
          T = 10. ** self.get('logT')[::-1]
          V = self.get('velocity')[::-1]
  
          # Abundances:
          def make_list(element,lowA,highA):
              l = []
              for i in range(lowA,highA+1):
                  l.append(element+str(i))
              return l
  
          abun_avail = list(self.cols.keys())
  
          def elemental_abund(ilist,abun_avail):
              X = np.zeros(len(self.get('mass')))
              for a in ilist:
                  if a in abun_avail:
                      X += self.get(a)[::-1]
  
              return X
  
          iH = ['h1','h2','prot']
          XH  = elemental_abund(iH, abun_avail)
          XHe = elemental_abund(make_list('he',1,5), abun_avail)
          XC  = elemental_abund(make_list('c',11,15), abun_avail)
          XN  = elemental_abund(make_list('n',12,16), abun_avail)
          XO  = elemental_abund(make_list('o',13,20), abun_avail)
          XNe = elemental_abund(make_list('ne',17,25), abun_avail)
          XNa = elemental_abund(make_list('na',20,25), abun_avail)
          XMg = elemental_abund(make_list('mg',21,28), abun_avail)
          XAl = elemental_abund(make_list('al',21,30), abun_avail)
          XSi = elemental_abund(make_list('si',25,34), abun_avail)
          XS  = elemental_abund(make_list('s',28,38), abun_avail)
          XAr = elemental_abund(make_list('ar',32,46), abun_avail)
          XCa = elemental_abund(make_list('ca',36,53), abun_avail)
          XFe = elemental_abund(make_list('fe',50,65), abun_avail)
          XCo = elemental_abund(make_list('co',52,66), abun_avail)
          XNi = elemental_abund(make_list('ni',54,71), abun_avail)
  
          XNi56 = self.get('ni56')
  
          # Write the output files:
          file_hyd = name+'.hyd'
          file_abn = name+'.abn'
  
          f = open(file_hyd,'w')
          # write header:
          f.write('  0.000E+00\n')
          f.write('# No.')
          f.write('Mr'.rjust(28)+
                  'dM'.rjust(28)+
                  'R'.rjust(28)+
                  'dR'.rjust(28)+
                  'Rho'.rjust(28)+
                  'PRE'.rjust(28)+
                  'T'.rjust(28)+
                  'V'.rjust(28)+
                  '\n')
          # write data:
          for i in range(len(zn)):
              f.write( str(zn[i]).rjust(5) +
                      '%.16E'.rjust(11) %Mr[i] +
                      '%.16E'.rjust(11) %dM[i] +
                      '%.16E'.rjust(11) %R[i] +
                      '%.16E'.rjust(11) %dR[i] +
                      '%.16E'.rjust(11) %Rho[i] +
                      '%.16E'.rjust(11) %PRE[i] +
                      '%.16E'.rjust(11) %T[i] +
                      '%.16E'.rjust(11) %V[i] +
                      '\n')
  
          f.close()
  
          f = open(file_abn,'w')
          # write header:
          f.write('# No.')
          f.write('Mr'.rjust(28)+
                  'H'.rjust(28)+
                  'He'.rjust(28)+
                  'C'.rjust(28)+
                  'N'.rjust(28)+
                  'O'.rjust(28)+
                  'Ne'.rjust(28)+
                  'Na'.rjust(28)+
                  'Mg'.rjust(28)+
                  'Al'.rjust(28)+
                  'Si'.rjust(28)+
                  'S'.rjust(28)+
                  'Ar'.rjust(28)+
                  'Ca'.rjust(28)+
                  'Fe'.rjust(28)+
                  'Co'.rjust(28)+
                  'Ni'.rjust(28)+
                  'X(56Ni)'.rjust(28)+
                  '\n')
          # write data:
          for i in range(len(zn)):
              f.write( str(zn[i]).rjust(5) +
                      '%.16E'.rjust(11) %Mr[i] +
                      '%.16E'.rjust(11) %XH[i] +
                      '%.16E'.rjust(11) %XHe[i] +
                      '%.16E'.rjust(11) %XC[i] +
                      '%.16E'.rjust(11) %XN[i] +
                      '%.16E'.rjust(11) %XO[i] +
                      '%.16E'.rjust(11) %XNe[i] +
                      '%.16E'.rjust(11) %XNa[i] +
                      '%.16E'.rjust(11) %XMg[i] +
                      '%.16E'.rjust(11) %XAl[i] +
                      '%.16E'.rjust(11) %XSi[i] +
                      '%.16E'.rjust(11) %XS[i] +
                      '%.16E'.rjust(11) %XAr[i] +
                      '%.16E'.rjust(11) %XCa[i] +
                      '%.16E'.rjust(11) %XFe[i] +
                      '%.16E'.rjust(11) %XCo[i] +
                      '%.16E'.rjust(11) %XNi[i] +
                      '%.16E'.rjust(11) %XNi56[i] +
                      '\n')