Index: trunk/SophyaLib/TArray/Makefile
===================================================================
--- trunk/SophyaLib/TArray/Makefile	(revision 3791)
+++ trunk/SophyaLib/TArray/Makefile	(revision 3809)
@@ -4,8 +4,8 @@
  
 clean:
-	 rm -f $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
+	 rm -f $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)fiospsqmtx.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)spesqmtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
 	rm -f $(SOPHYALIBP)libTArray.a
  
-$(SOPHYALIBP)libTArray.a : $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
+$(SOPHYALIBP)libTArray.a : $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)fiospsqmtx.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)spesqmtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
 	$(AR) $(ARFLAGS) $@ $? 
 	touch $(SOPHYAINCP)/SophyaConfInfo/libsophya.objlist
@@ -13,19 +13,25 @@
 $(SOPHYAOBJP)basarr.o: basarr.cc $(SOPHYAINCP)sopnamsp.h \
   $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
   $(SOPHYAINCP)pexceptions.h \
   $(SOPHYAINCP)machdefs.h basarr.h \
   $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h
 	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  basarr.cc
  
@@ -33,171 +39,309 @@
   $(SOPHYAINCP)pexceptions.h \
   $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
   $(SOPHYAINCP)fiondblock.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)ndatablock.h \
-  $(SOPHYAINCP)anydataobj.h fioarr.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h fioarr.h tarray.h basarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  utilarr.h $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tmatrix_tsnl.h \
+  tvector.h tvector_tsnl.h $(SOPHYAINCP)datatype.h \
+  $(SOPHYAINCP)machdefs.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  fioarr.cc
+ 
+$(SOPHYAOBJP)fiospsqmtx.o: fiospsqmtx.cc $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)fiondblock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h fiospsqmtx.h spesqmtx.h utilarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tarray.h basarr.h \
+  $(SOPHYAINCP)dvlist.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)mutyv.h tmatrix_tsnl.h trngmtx.h \
+  spesqmtx_tsnl.h diagmtx.h symmtx.h $(SOPHYAINCP)datatype.h \
+  $(SOPHYAINCP)machdefs.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  fiospsqmtx.cc
+ 
+$(SOPHYAOBJP)matharr.o: matharr.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h matharr.h tarray.h basarr.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  matharr.cc
+ 
+$(SOPHYAOBJP)sopemtx.o: sopemtx.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h sopemtx.h tmatrix.h tarray.h \
+  basarr.h $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix_tsnl.h tvector.h \
+  tvector_tsnl.h $(SOPHYAINCP)smathconst.h \
+  $(SOPHYAINCP)machdefs.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  sopemtx.cc
+ 
+$(SOPHYAOBJP)spesqmtx.o: spesqmtx.cc spesqmtx.h $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h utilarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tarray.h basarr.h \
+  $(SOPHYAINCP)dvlist.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)mutyv.h tmatrix_tsnl.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  spesqmtx.cc
+ 
+$(SOPHYAOBJP)tarray.o: tarray.cc $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)pexceptions.h \
   $(SOPHYAINCP)machdefs.h tarray.h basarr.h \
   $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix.h \
-  tmatrix_tsnl.h tvector.h tvector_tsnl.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)datatype.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  fioarr.cc
- 
-$(SOPHYAOBJP)matharr.o: matharr.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h matharr.h tarray.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarray.cc
+ 
+$(SOPHYAOBJP)tarrinit.o: tarrinit.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h tarrinit.h \
+  $(SOPHYAINCP)sophyainit.h \
+  $(SOPHYAINCP)machdefs.h fioarr.h tarray.h basarr.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tmatrix_tsnl.h \
+  tvector.h tvector_tsnl.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarrinit.cc
+ 
+$(SOPHYAOBJP)tmatrix.o: tmatrix.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h tmatrix.h tarray.h basarr.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix_tsnl.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tmatrix.cc
+ 
+$(SOPHYAOBJP)tvector.o: tvector.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h tvector.h tmatrix.h tarray.h \
   basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  matharr.cc
- 
-$(SOPHYAOBJP)sopemtx.o: sopemtx.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h sopemtx.h tmatrix.h \
-  tarray.h basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix_tsnl.h \
-  tvector.h tvector_tsnl.h \
-  $(SOPHYAINCP)smathconst.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  sopemtx.cc
- 
-$(SOPHYAOBJP)tarray.o: tarray.cc $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)machdefs.h tarray.h basarr.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarray.cc
- 
-$(SOPHYAOBJP)tarrinit.o: tarrinit.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h tarrinit.h \
-  $(SOPHYAINCP)sophyainit.h \
-  $(SOPHYAINCP)machdefs.h fioarr.h tarray.h \
-  basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix.h \
-  tmatrix_tsnl.h tvector.h tvector_tsnl.h \
-  $(SOPHYAINCP)ppersist.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarrinit.cc
- 
-$(SOPHYAOBJP)tmatrix.o: tmatrix.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)machdefs.h tmatrix.h tarray.h \
-  basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix_tsnl.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tmatrix.cc
- 
-$(SOPHYAOBJP)tvector.o: tvector.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)machdefs.h tvector.h tmatrix.h \
-  tarray.h basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix_tsnl.h \
-  tvector_tsnl.h
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix_tsnl.h tvector_tsnl.h
 	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tvector.cc
  
 $(SOPHYAOBJP)utilarr.o: utilarr.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h utilarr.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h utilarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
   $(SOPHYAINCP)randr48.h \
-  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
   $(SOPHYAINCP)randinterf.h
 	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  utilarr.cc
Index: trunk/SophyaLib/TArray/array.h
===================================================================
--- trunk/SophyaLib/TArray/array.h	(revision 3791)
+++ trunk/SophyaLib/TArray/array.h	(revision 3809)
@@ -1,6 +1,7 @@
-// This may look like C code, but it is really -*- C++ -*-
+// This may look like C code, but it is really -*- C++ -*-    
+//   --- SOPHYA software - TArray module ---
+//   (C) UPS+LAL IN2P3/CNRS     (C) DAPNIA-SPP/CEA 
 //  All the necessari include files for array manipulation
-//                     R. Ansari, C.Magneville   03/2000
-
+                    
 #ifndef AllTArray_SEEN
 #define AllTArray_SEEN
@@ -13,4 +14,8 @@
 #include "sopemtx.h"
 #include "arrctcast.h"
+#include "spesqmtx.h"
+#include "diagmtx.h"
+#include "trngmtx.h"
+#include "symmtx.h"
 
 #endif
Index: trunk/SophyaLib/TArray/diagmtx.h
===================================================================
--- trunk/SophyaLib/TArray/diagmtx.h	(revision 3809)
+++ trunk/SophyaLib/TArray/diagmtx.h	(revision 3809)
@@ -0,0 +1,211 @@
+// This may look like C code, but it is really -*- C++ -*-
+
+#ifndef DIAGMTX_H_SEEN
+#define DIAGMTX_H_SEEN
+
+#include "spesqmtx.h"
+
+// doit etre mis en dehors du namespace
+/*!
+  \class SOPHYA::DiagonalMatrix
+  \ingroup TArray
+  \brief Class representing a diagonal matrix. 
+
+  This class offers similar functionalities to the TArray<T> / TMatrix<T> classes, like 
+  reference sharing and counting, arithmetic operators ... However, this class has no 
+  sub matrix extraction method.
+*/
+
+namespace SOPHYA {
+
+//! Class for inferior triangular matrix (base class for the class Alm)
+template <class T>
+class DiagonalMatrix : public SpecialSquareMatrix<T> {
+public :
+
+#include "spesqmtx_tsnl.h"
+
+//! Default constructor - TriangMatrix of size 0, SetSize() should be called before the object is used
+explicit DiagonalMatrix()
+  : SpecialSquareMatrix<T>(C_DiagonalMatrix) 
+{ 
+  mOffDiag = T(0);
+}
+
+//! Instanciate a triangular matrix from the number of rows (rowSize must be > 0)
+explicit DiagonalMatrix(sa_size_t rowSize)
+  : SpecialSquareMatrix<T>(rowSize, C_DiagonalMatrix)
+{
+  if (rowSize < 1) 
+    throw ParmError("DiagonalMatrix<T>::DiagonalMatrix(rsz) rsz <= 0");
+  mElems.ReSize(rowSize); 
+  mInfo = NULL;
+  mOffDiag = T(0);
+}
+
+//! Copy constructor (possibility of sharing datas)
+DiagonalMatrix(DiagonalMatrix<T> const & a,  bool share=false)  
+  : SpecialSquareMatrix<T>(a, share) 
+{
+  mOffDiag = T(0);
+}
+
+//! Copy constructor (possibility of sharing datas)
+DiagonalMatrix(SpecialSquareMatrix<T> const & a,  bool share=false)  
+  : SpecialSquareMatrix<T>(a, share) 
+{
+  if (a.MtxType() != C_DiagonalMatrix) 
+    throw TypeMismatchExc("DiagonalMatrix(a) a NOT a DiagonalMatrix");
+  mOffDiag = T(0);
+}
+
+/*!
+ \brief Create a lower triangular matrix from a square matrix. 
+ Off diagonal elements are ignored.
+*/
+explicit DiagonalMatrix(TMatrix<T> const & mx)  
+  : SpecialSquareMatrix<T>(C_DiagonalMatrix)
+{
+  if ((mx.NRows() != mx.NCols()) || (mx.NRows() < 1)) 
+    throw ParmError("DiagonalMatrix<T>::(TMatrix<T> const & mx) mx not allocated OR NOT a square matrix");
+  SetSize(mx.NRows());
+  for(sa_size_t l=0; l<NRows(); l++) (*this)(l,l) = mx(l,l);
+  mOffDiag = T(0);
+}
+
+//! Sets or change the triangular matrix size, specifying the new number of rows
+virtual void SetSize(sa_size_t rowSize) 
+{
+  if (rowSize < 1) 
+    throw ParmError("DiagonalMatrix<T>::SetSize(rsz) rsz <= 0");
+  if (rowSize == mNrows)  return;
+  mNrows=rowSize;
+  mElems.ReSize(mNrows);  
+}
+
+
+//! Return the object (diagonal matrix) as a standard (TMatrix<T>) square matrix
+virtual TMatrix<T> ConvertToStdMatrix() const
+{
+  if (mNrows < 1)
+    throw SzMismatchError("DiagonalMatrix<T>::ConvertToStdMatrix() (this) not allocated !");
+  TMatrix<T> mx(NRows(), NRows());
+  for(sa_size_t l=0; l<NRows(); l++) mx(l,l) = (*this)(l,l);
+  return mx;
+}
+ 
+
+//--- Operateurs = (T b) , = (DiagonalMatrix<T> b)
+//! operator = a , to set all elements to the value \b a
+inline DiagonalMatrix<T>& operator = (T a)
+  {  SetCst(a);  return (*this);  }
+//! operator = DiagonalMatrix<T> a , element by element copy operator
+inline DiagonalMatrix<T>& operator = (DiagonalMatrix<T> const & a)
+  {  Set(a); return (*this); }
+//! operator = Sequence seq  
+inline DiagonalMatrix<T>& operator = (Sequence const & seq)
+  {  SetSeq(seq); return (*this); }
+//! operator = Sequence seq  
+inline DiagonalMatrix<T>& operator = (IdentityMatrix & idmx)
+  {  SetCst(idmx.Diag()); return (*this); }
+
+//--- Operateurs d'acces aux elements 
+//! Element access operator (R/W): access to elements row \b r and column \b c
+inline T& operator()(sa_size_t r, sa_size_t c) 
+{ 
+  if ((r<0)||(r>=mNrows))  
+    throw RangeCheckError("DiagonalMatrix<T>::operator()(r,c) (r<0)||(r>=NRows())");
+  if (r!=c) { mOffDiag = T(0); return mOffDiag; }
+  return mElems(r);
+}
+//! Element access operator (RO): access to elements row \b r and column \b c
+inline T operator()(sa_size_t r, sa_size_t c) const
+{
+  if ((r<0)||(r>=mNrows))  
+    throw RangeCheckError("DiagonalMatrix<T>::operator()(r,c) (r<0)||(r>=NRows())");
+  if (r!=c) { mOffDiag = T(0); return mOffDiag; }
+  return mElems(r);
+}
+
+//! Diagonal Matrix product (multiplication) : ret_matrix = (*this) * dmx
+DiagonalMatrix<T> Multiply(DiagonalMatrix<T> const & dmx) const
+{
+  if (NRows() != dmx.NRows()) 
+    throw SzMismatchError("DiagonalMatrix<T>::Multiply(DiagonalMatrix<T> dmx): different sizes");
+  DiagonalMatrix<T> ret(NRows());
+  for(size_t k=0; k<mElems.Size(); k++) 
+    ret.mElems(k) = mElems(k)*dmx.mElems(k);
+  ret.SetTemp(true);
+  return ret;
+}
+
+//! Matrix product (multiplication) : ret_matrix = (*this) * mx
+TMatrix<T> MultiplyDG(TMatrix<T> const & mx) const 
+{
+  if (NCols() != mx.NRows())
+    throw SzMismatchError("DiagonalMatrix<T>::MultiplyDG(TMatrix<T> mx): NCols()!=mx.NRows()");
+    
+  TMatrix<T> ret(mx, false);
+  for(sa_size_t r=0; r<NRows(); r++) 
+    ret.Row(r) *= (*this)(r,r);
+  ret.SetTemp(true);
+  return ret;
+}
+
+//! Matrix product (multiplication) : ret_matrix = mx * (*this)
+TMatrix<T> MultiplyGD(TMatrix<T> const & mx) const
+{
+  if (NRows() != mx.NCols())
+    throw SzMismatchError("DiagonalMatrix<T>::MultiplyGD(TMatrix<T> mx): NRows()!=mx.NCols()");
+    
+  TMatrix<T> ret(mx, false);
+  for(sa_size_t c=0; c<NRows(); c++) 
+    ret.Column(c) *= (*this)(c,c);
+  ret.SetTemp(true);
+  return ret;
+}
+
+
+//! ASCII dump/print of the triangular matrix object (nprt=-1 for printing all diagonal elements)
+ostream& Print(ostream& os, sa_size_t nprt=-1) const
+{
+  os << "DiagonalMatrix< " << typeid(T).name() 
+     << " > NRow=" << mNrows << " NbElem<>0 : " << Size() << endl;
+  for(sa_size_t k=0; k<Size(); k+=8) {
+    if (k%32==0) os << "DiagonalElements: [ " << k << " ..." << k+31 <<" ] :" << endl;
+    sa_size_t jmx=k+8;
+    if (jmx>Size()) jmx = Size();
+    for(sa_size_t j=k; j<jmx; j++)  os << mElems(j) << " , ";
+    os << endl; 
+    if (k >= nprt)  break;
+  } 
+  return os;
+}
+
+protected: 
+mutable T mOffDiag;
+};
+
+//----- Surcharge d'operateurs C = A * B (multiplication matricielle)
+/*! \ingroup TArray \fn operator*(const DiagonalMatrix<T>&,const DiagonalMatrix<T>&)
+  \brief * : DiagonalMatrix multiplication \b a and \b b */
+template <class T> 
+inline DiagonalMatrix<T> operator * (const DiagonalMatrix<T>& a, const DiagonalMatrix<T>& b)
+   { return(a.Multiply(b)); }
+
+/*! \ingroup TArray \fn operator*(const DiagonalMatrix<T>&,const TMatrix<T>&)
+  \brief * : Matrix multiplication DiagonalMatrix (\b a ) *  TMatrix<T> ( \b b ) */
+template <class T> 
+inline TMatrix<T> operator * (const DiagonalMatrix<T>& a, const TMatrix<T>& b)
+   { return(a.MultiplyDG(b)); }
+
+/*! \ingroup TArray \fn operator*(const TMatrix<T>&,const DiagonalMatrix<T>&)
+  \brief * : Matrix multiplication TMatrix (\b a ) *  DiagonalMatrix<T> ( \b b ) */
+template <class T> 
+inline TMatrix<T> operator * (const TMatrix<T>& a, const DiagonalMatrix<T>& b)
+   { return(b.MultiplyGD(a)); }
+
+  
+}   // namespace SOPHYA
+
+#endif
Index: trunk/SophyaLib/TArray/fiospsqmtx.cc
===================================================================
--- trunk/SophyaLib/TArray/fiospsqmtx.cc	(revision 3809)
+++ trunk/SophyaLib/TArray/fiospsqmtx.cc	(revision 3809)
@@ -0,0 +1,204 @@
+//     Persistence manager for template numerical arrays
+//               R. Ansari, C.Magneville   03/2000
+
+#include "pexceptions.h"
+#include "fiondblock.h"
+#include "fiospsqmtx.h"
+#include "datatype.h"
+#include <typeinfo>
+
+
+//----------------------------------------------------------------------
+/*!
+  \class SOPHYA::FIO_SpecialSquareMatrix
+  \ingroup TArray
+
+  \brief Manager class the PPF I/O (SOPHYA persistence) of SpecialSquareMatrix<T> objects.
+ 
+ */
+///////////////////////////////////////////////////////////
+
+namespace SOPHYA {
+
+//! Default constructor
+template <class T>
+FIO_SpecialSquareMatrix<T>::FIO_SpecialSquareMatrix()
+{
+  dobj=NULL;
+  ownobj=false;
+}
+
+
+//! Constructor from the file \b filename - Reads the SpecialSquareMatrix object from the specified file.
+template <class T>
+FIO_SpecialSquareMatrix<T>::FIO_SpecialSquareMatrix(string const & filename) 
+{
+  dobj=NULL;
+  ownobj=false;
+  Read(filename);
+}
+
+//! Constructor from the SpecialSquareMatrix \b obj
+template <class T>
+FIO_SpecialSquareMatrix<T>::FIO_SpecialSquareMatrix(const SpecialSquareMatrix<T> & obj) 
+{ 
+  const SymmetricMatrix<T>* sym = dynamic_cast< const SymmetricMatrix<T> * >(&obj);
+  if (sym != NULL) dobj = new SymmetricMatrix<T>(*sym, true); 
+  else {
+    const LowerTriangularMatrix<T>* trn = dynamic_cast< const LowerTriangularMatrix<T> * >(&obj);
+    if (trn != NULL)  dobj = new LowerTriangularMatrix<T>(*trn, true);
+    else {
+      const DiagonalMatrix<T>* diag = dynamic_cast< const DiagonalMatrix<T> * >(&obj);
+      if (diag != NULL)  dobj = new DiagonalMatrix<T>(*diag, true);
+      else throw TypeMismatchExc("FIO_SpecialSquareMatrix<T>::FIO_SpecialSquareMatrix(obj) Not a TriangularMatrix<T> / DiagonalMatrix<T> / SymmetricMatrix<T> ");
+    }
+  }
+  ownobj=true; 
+}
+
+//! Connect with a SpecialSquareMatrix \b obj
+template <class T>
+FIO_SpecialSquareMatrix<T>::FIO_SpecialSquareMatrix(SpecialSquareMatrix<T> * obj) 
+{ 
+  dobj = obj;
+  ownobj=false; 
+}
+
+//! destructor
+template <class T>
+FIO_SpecialSquareMatrix<T>::~FIO_SpecialSquareMatrix()
+{
+  if (ownobj && dobj) delete dobj;
+}
+
+//! Return the pointer to the connected SpecialSquareMatrix<T> object
+template <class T>
+AnyDataObj* FIO_SpecialSquareMatrix<T>::DataObj()
+{
+  return(dobj);
+}
+
+//! Connect SpecialSquareMatrix \b o
+template <class T>
+void FIO_SpecialSquareMatrix<T>::SetDataObj(AnyDataObj & o)
+{
+  SpecialSquareMatrix<T> * po = dynamic_cast< SpecialSquareMatrix<T> * >(&o);
+  if (po == NULL) {
+    char buff[160];
+    sprintf(buff,"FIO_SpecialSquareMatrix<%s>::SetDataObj(%s) - Object type  error ! ",
+	    DataTypeInfo<T>::getTypeName().c_str(), typeid(o).name());
+    throw TypeMismatchExc(PExcLongMessage(buff));    
+  }
+  if (ownobj && dobj) delete dobj;
+  dobj = po; ownobj = false;
+} 
+
+//! Reads in the SpecialSquareMatrix<T> object from the PPF input stream \b is
+template <class T>
+void FIO_SpecialSquareMatrix<T>::ReadSelf(PInPersist& is)
+{
+// On lit les 4 premiers uint_4
+//  0: Numero de version,  1 : MtxType,  2 != 0 , has Info, 3 : Reserved
+  uint_4 itab[4];
+  is.Get(itab,4);
+
+  if (dobj == NULL)  { 
+    if (itab[1] == (uint_4)C_DiagonalMatrix)  dobj = new DiagonalMatrix<T>;
+    else if (itab[1] == (uint_4)C_LowerTriangularMatrix)  dobj = new LowerTriangularMatrix<T>;
+    else if (itab[1] == (uint_4)C_SymmetricMatrix)  dobj = new SymmetricMatrix<T>;
+    else throw TypeMismatchExc("FIO_SpecialSquareMatrix<T>::ReadSelf() bad matrix type in file");
+  }
+  else {
+    if ((uint_4)dobj->MtxType() != itab[1]) {
+      // cout << " BUG*** dobj->MtxType()=" << (uint_4)(dobj->MtxType()) << " itab[1]=" << itab[1] << endl;
+      throw TypeMismatchExc("FIO_SpecialSquareMatrix<T>::ReadSelf() matrix type dobj <> in file ");
+    }
+  }
+
+// On lit les tailles, etc ...
+
+// Tous les sa_size_t sont ecrit/lu en int_8 afin de maintenir la compatibilite  
+// -- 0: rowSize
+  int_8 tmpi8s[3];
+
+  is.Get(tmpi8s, 3);
+  dobj->SetSize(tmpi8s[0]);
+
+// On lit le datablock
+  is >> dobj->DataBlock();
+// On lit le DVList info si necessaire
+  if (itab[2] != 0)  is >> dobj->Info();  
+}
+
+//! Writes the SpecialSquareMatrix<T> object to the output PPF  stream \b os
+template <class T>
+void FIO_SpecialSquareMatrix<T>::WriteSelf(POutPersist& os) const
+{
+  if (dobj == NULL)   return;
+//  On ecrit 3 uint_4 .... 
+//  0: Numero de version,  1 : reserved  2 != 0 , has Info
+  uint_4 itab[4];
+  itab[0] = 1;  // Numero de version a 1 - fev 2009
+  itab[1] = (uint_4)dobj->MtxType();  // Type de matrice : Diagonal, Triangular, Symmetric ...
+  // cout << " *** DBG** WriteSelf() dobj->MtxType()= " << itab[1] << endl;
+  itab[2] = (dobj->mInfo != NULL) ? 1 : 0;
+  itab[3] = 0;  // Reserved
+  os.Put(itab,4);
+
+// On ecrit les tailles, etc ...
+// Tous les sa_size_t sont ecrit en int_8 afin de pouvoir etre ecrit/relu
+// entre version du programme utilisant int_4 OU int_8 pour sa_size_t 
+  int_8 tmpi8s[3];
+  tmpi8s[0] = dobj->NRows();
+  tmpi8s[1] = tmpi8s[2] = 0;
+  os.Put(tmpi8s, 3);
+  // On ecrit le datablock
+  os << dobj->DataBlock();
+// On ecrit le DVList info si necessaire
+  if (itab[2] != 0)  os << dobj->Info();
+}
+
+} // Fin de namespace SOPHYA 
+
+
+///////////////////////////////////////////////////////////////
+#ifdef __CXX_PRAGMA_TEMPLATES__
+// Instances des delegues FileIO (PPersist)
+#pragma define_template FIO_SpecialSquareMatrix<uint_1>
+#pragma define_template FIO_SpecialSquareMatrix<uint_2>
+#pragma define_template FIO_SpecialSquareMatrix<uint_4>
+#pragma define_template FIO_SpecialSquareMatrix<uint_8>
+#pragma define_template FIO_SpecialSquareMatrix<int_1>
+#pragma define_template FIO_SpecialSquareMatrix<int_2>
+#pragma define_template FIO_SpecialSquareMatrix<int_4>
+#pragma define_template FIO_SpecialSquareMatrix<int_8>
+#pragma define_template FIO_SpecialSquareMatrix<r_4>
+#pragma define_template FIO_SpecialSquareMatrix<r_8>
+#pragma define_template FIO_SpecialSquareMatrix< complex<r_4> >
+#pragma define_template FIO_SpecialSquareMatrix< complex<r_8> >
+#ifdef SO_LDBLE128
+#pragma define_template FIO_SpecialSquareMatrix<r_16>
+#pragma define_template FIO_SpecialSquareMatrix< complex<r_16> >
+#endif
+#endif
+
+#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
+namespace SOPHYA {
+template class FIO_SpecialSquareMatrix<uint_1>;
+template class FIO_SpecialSquareMatrix<uint_2>;
+template class FIO_SpecialSquareMatrix<uint_4>;
+template class FIO_SpecialSquareMatrix<uint_8>;
+template class FIO_SpecialSquareMatrix<int_1>;
+template class FIO_SpecialSquareMatrix<int_2>;
+template class FIO_SpecialSquareMatrix<int_4>;
+template class FIO_SpecialSquareMatrix<int_8>;
+template class FIO_SpecialSquareMatrix<r_4>;
+template class FIO_SpecialSquareMatrix<r_8>;
+template class FIO_SpecialSquareMatrix< complex<r_4> >;
+template class FIO_SpecialSquareMatrix< complex<r_8> >;
+#ifdef SO_LDBLE128
+template class FIO_SpecialSquareMatrix<r_16>;
+template class FIO_SpecialSquareMatrix< complex<r_16> >;
+#endif
+}
+#endif
Index: trunk/SophyaLib/TArray/fiospsqmtx.h
===================================================================
--- trunk/SophyaLib/TArray/fiospsqmtx.h	(revision 3809)
+++ trunk/SophyaLib/TArray/fiospsqmtx.h	(revision 3809)
@@ -0,0 +1,61 @@
+// This may look like C code, but it is really -*- C++ -*-
+//      template array class for numerical types
+//                     R. Ansari, C.Magneville   03/2000
+//  Ajout gestionnaire PPF pour TriangularMatrix<T>   02/2009
+
+#ifndef FIOSpSqMtx_SEEN
+#define FIOSpSqMtx_SEEN
+
+#include "machdefs.h"
+#include <iostream>
+#include "spesqmtx.h"
+#include "trngmtx.h"
+#include "diagmtx.h"
+#include "symmtx.h"
+#include "ppersist.h"
+
+namespace SOPHYA {
+
+//------------------------------------------------------------------------------------
+
+/////////////////////////////////////////////////////////////////////////
+// Gestionnaire PPF de SpecialSquareMatrix<T>
+
+/*! 
+  \brief Manager class the PPF I/O (SOPHYA persistence) of TriangularMatrix<T> objects.
+ */
+
+template <class T>
+class FIO_SpecialSquareMatrix : public  PPersist  {
+public:
+  FIO_SpecialSquareMatrix();
+  FIO_SpecialSquareMatrix(string const & filename); 
+  FIO_SpecialSquareMatrix(const SpecialSquareMatrix<T> & obj);
+  FIO_SpecialSquareMatrix(SpecialSquareMatrix<T> * obj);
+  virtual ~FIO_SpecialSquareMatrix();
+  virtual AnyDataObj* DataObj();
+  virtual void        SetDataObj(AnyDataObj & o);
+  inline operator SpecialSquareMatrix<T>() { return(*dobj); }
+protected :
+  virtual void ReadSelf(PInPersist&);           
+  virtual void WriteSelf(POutPersist&) const;  
+  SpecialSquareMatrix<T> * dobj; //!< Pointer to the connected SpecialSquareMatrix
+  bool ownobj; //!< true if the connected object belong to the class
+};
+
+/*! \ingroup TArray \fn operator<<(POutPersist&,SpecialSquareMatrix<T>&)
+  \brief Write SpecialSquareMatrix \b obj into POutPersist stream \b os */
+template <class T>
+inline POutPersist& operator << (POutPersist& os, SpecialSquareMatrix<T> & obj)
+{ FIO_SpecialSquareMatrix<T> fio(&obj);  fio.Write(os);  return(os); }
+
+/*! \ingroup TArray \fn operator>>(PInPersist&,SpecialSquareMatrix<T>&)
+  \brief Read SpecialSquareMatrix \b obj from PInPersist stream \b os */
+template <class T>
+inline PInPersist& operator >> (PInPersist& is, SpecialSquareMatrix<T> & obj)
+{ FIO_SpecialSquareMatrix<T> fio(&obj); is.SkipToNextObject(); fio.Read(is); return(is); }
+
+
+} // Fin du namespace
+
+#endif
Index: trunk/SophyaLib/TArray/objlist.list
===================================================================
--- trunk/SophyaLib/TArray/objlist.list	(revision 3791)
+++ trunk/SophyaLib/TArray/objlist.list	(revision 3809)
@@ -1,6 +1,8 @@
 basarr.o
 fioarr.o
+fiospsqmtx.o
 matharr.o
 sopemtx.o
+spesqmtx.o
 tarray.o
 tarrinit.o
Index: trunk/SophyaLib/TArray/smakefile
===================================================================
--- trunk/SophyaLib/TArray/smakefile	(revision 3791)
+++ trunk/SophyaLib/TArray/smakefile	(revision 3809)
@@ -4,8 +4,8 @@
  
 clean:
-	 rm -f $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
+	 rm -f $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)fiospsqmtx.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)spesqmtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
 	rm -f $(SOPHYALIBP)libTArray.a
  
-$(SOPHYALIBP)libTArray.a : $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
+$(SOPHYALIBP)libTArray.a : $(SOPHYAOBJP)basarr.o $(SOPHYAOBJP)fioarr.o $(SOPHYAOBJP)fiospsqmtx.o $(SOPHYAOBJP)matharr.o $(SOPHYAOBJP)sopemtx.o $(SOPHYAOBJP)spesqmtx.o $(SOPHYAOBJP)tarray.o $(SOPHYAOBJP)tarrinit.o $(SOPHYAOBJP)tmatrix.o $(SOPHYAOBJP)tvector.o $(SOPHYAOBJP)utilarr.o
 	$(AR) $(ARFLAGS) $@ $? 
 	touch $(SOPHYAINCP)/SophyaConfInfo/libsophya.objlist
@@ -13,19 +13,25 @@
 $(SOPHYAOBJP)basarr.o: basarr.cc $(SOPHYAINCP)sopnamsp.h \
   $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
   $(SOPHYAINCP)pexceptions.h \
   $(SOPHYAINCP)machdefs.h basarr.h \
   $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h
 	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  basarr.cc
  
@@ -33,171 +39,309 @@
   $(SOPHYAINCP)pexceptions.h \
   $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
   $(SOPHYAINCP)fiondblock.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)ndatablock.h \
-  $(SOPHYAINCP)anydataobj.h fioarr.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h fioarr.h tarray.h basarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  utilarr.h $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tmatrix_tsnl.h \
+  tvector.h tvector_tsnl.h $(SOPHYAINCP)datatype.h \
+  $(SOPHYAINCP)machdefs.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  fioarr.cc
+ 
+$(SOPHYAOBJP)fiospsqmtx.o: fiospsqmtx.cc $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)fiondblock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h fiospsqmtx.h spesqmtx.h utilarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tarray.h basarr.h \
+  $(SOPHYAINCP)dvlist.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)mutyv.h tmatrix_tsnl.h trngmtx.h \
+  spesqmtx_tsnl.h diagmtx.h symmtx.h $(SOPHYAINCP)datatype.h \
+  $(SOPHYAINCP)machdefs.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  fiospsqmtx.cc
+ 
+$(SOPHYAOBJP)matharr.o: matharr.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h matharr.h tarray.h basarr.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  matharr.cc
+ 
+$(SOPHYAOBJP)sopemtx.o: sopemtx.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h sopemtx.h tmatrix.h tarray.h \
+  basarr.h $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix_tsnl.h tvector.h \
+  tvector_tsnl.h $(SOPHYAINCP)smathconst.h \
+  $(SOPHYAINCP)machdefs.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  sopemtx.cc
+ 
+$(SOPHYAOBJP)spesqmtx.o: spesqmtx.cc spesqmtx.h $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h utilarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tarray.h basarr.h \
+  $(SOPHYAINCP)dvlist.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)mutyv.h tmatrix_tsnl.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  spesqmtx.cc
+ 
+$(SOPHYAOBJP)tarray.o: tarray.cc $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)pexceptions.h \
   $(SOPHYAINCP)machdefs.h tarray.h basarr.h \
   $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix.h \
-  tmatrix_tsnl.h tvector.h tvector_tsnl.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)datatype.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  fioarr.cc
- 
-$(SOPHYAOBJP)matharr.o: matharr.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h matharr.h tarray.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarray.cc
+ 
+$(SOPHYAOBJP)tarrinit.o: tarrinit.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h tarrinit.h \
+  $(SOPHYAINCP)sophyainit.h \
+  $(SOPHYAINCP)machdefs.h fioarr.h tarray.h basarr.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix.h tmatrix_tsnl.h \
+  tvector.h tvector_tsnl.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarrinit.cc
+ 
+$(SOPHYAOBJP)tmatrix.o: tmatrix.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h tmatrix.h tarray.h basarr.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix_tsnl.h
+	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tmatrix.cc
+ 
+$(SOPHYAOBJP)tvector.o: tvector.cc $(SOPHYAINCP)sopnamsp.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h tvector.h tmatrix.h tarray.h \
   basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  matharr.cc
- 
-$(SOPHYAOBJP)sopemtx.o: sopemtx.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h sopemtx.h tmatrix.h \
-  tarray.h basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix_tsnl.h \
-  tvector.h tvector_tsnl.h \
-  $(SOPHYAINCP)smathconst.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  sopemtx.cc
- 
-$(SOPHYAOBJP)tarray.o: tarray.cc $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)machdefs.h tarray.h basarr.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarray.cc
- 
-$(SOPHYAOBJP)tarrinit.o: tarrinit.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h tarrinit.h \
-  $(SOPHYAINCP)sophyainit.h \
-  $(SOPHYAINCP)machdefs.h fioarr.h tarray.h \
-  basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix.h \
-  tmatrix_tsnl.h tvector.h tvector_tsnl.h \
-  $(SOPHYAINCP)ppersist.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tarrinit.cc
- 
-$(SOPHYAOBJP)tmatrix.o: tmatrix.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)machdefs.h tmatrix.h tarray.h \
-  basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix_tsnl.h
-	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tmatrix.cc
- 
-$(SOPHYAOBJP)tvector.o: tvector.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)machdefs.h tvector.h tmatrix.h \
-  tarray.h basarr.h $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)dvlist.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)ndatablock.h utilarr.h \
-  $(SOPHYAINCP)randinterf.h tmatrix_tsnl.h \
-  tvector_tsnl.h
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h $(SOPHYAINCP)dvlist.h \
+  $(SOPHYAINCP)objfio.h $(SOPHYAINCP)mutyv.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h utilarr.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h tmatrix_tsnl.h tvector_tsnl.h
 	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  tvector.cc
  
 $(SOPHYAOBJP)utilarr.o: utilarr.cc $(SOPHYAINCP)sopnamsp.h \
-  $(SOPHYAINCP)machdefs.h utilarr.h \
-  $(SOPHYAINCP)mutyv.h \
-  $(SOPHYAINCP)machdefs.h \
-  $(SOPHYAINCP)timestamp.h \
-  $(SOPHYAINCP)objfio.h \
-  $(SOPHYAINCP)anydataobj.h \
-  $(SOPHYAINCP)ppersist.h \
-  $(SOPHYAINCP)pexceptions.h \
-  $(SOPHYAINCP)gnumd5.h \
-  $(SOPHYAINCP)ppfbinstream.h \
-  $(SOPHYAINCP)rawstream.h \
-  $(SOPHYAINCP)ppfnametag.h \
-  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)osx_values.h utilarr.h \
+  $(SOPHYAINCP)mutyv.h $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)timestamp.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)pexceptions.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)gnumd5.h \
+  $(SOPHYAINCP)ppfbinstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)rawstream.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)ppfnametag.h \
+  $(SOPHYAINCP)ppersist.h \
+  $(SOPHYAINCP)randinterf.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
   $(SOPHYAINCP)randr48.h \
-  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h $(SOPHYAINCP)objfio.h \
+  $(SOPHYAINCP)ndatablock.h \
+  $(SOPHYAINCP)machdefs.h \
+  $(SOPHYAINCP)anydataobj.h \
   $(SOPHYAINCP)randinterf.h
 	$(CXXCOMPILE) $(CXXTEMPFLG) -o $@  utilarr.cc
Index: trunk/SophyaLib/TArray/spesqmtx.cc
===================================================================
--- trunk/SophyaLib/TArray/spesqmtx.cc	(revision 3809)
+++ trunk/SophyaLib/TArray/spesqmtx.cc	(revision 3809)
@@ -0,0 +1,316 @@
+// --- Special Square matrices base class 
+// This code is part of the SOPHYA library
+//  (C) Univ. Paris-Sud   (C) LAL-IN2P3/CNRS   (C) IRFU-CEA
+//  (C) R. Ansari, C.Magneville    2009  
+
+#include "spesqmtx.h"
+#include "pexceptions.h"
+#include "tmatrix.h"
+
+// doit etre mis en dehors du namespace
+/*!
+  \class SOPHYA::SpecialSquareMatrix
+  \ingroup TArray
+  \brief The base class for representing and manipulating special square matrices
+  Diagonal matrix, Symmetric matrix, Triangular matrix
+
+  This class implements the common operations for special square matrices. Only objects 
+  from the sub-classes (DiagonalMatrix<T>, LowerTriangularMatrix<T>, SymmetricMatrix<T>) 
+  can be instanciated. Theses classes offer similar functionalities to the TArray<T> / TMatrix<T> 
+  classes, like reference sharing and arithmetic operations. However, no sub-matrix extraction 
+  method is provided. 
+  sub matrix extraction method.
+*/
+
+namespace SOPHYA {
+
+//! Default constructor - SpecialSquareMatrix of size 0, SetSize() should be called before the object is used
+template <class T>
+SpecialSquareMatrix<T>::SpecialSquareMatrix(SpSqMtxType typ)
+  : mType(typ) , mNrows(0) , mInfo(NULL) 
+{
+}
+
+//! Instanciate a  SpecialSquareMatrix matrix from the number of rows (rowSize must be > 0)
+template <class T>
+SpecialSquareMatrix<T>::SpecialSquareMatrix(sa_size_t rowSize, SpSqMtxType typ)
+  : mType(typ) , mNrows(rowSize) , mInfo(NULL)
+{
+  if (rowSize < 1) 
+    throw ParmError("SpecialSquareMatrix<T>::SpecialSquareMatrix(rsz) rsz <= 0");
+}
+
+//! Copy constructor (possibility of sharing datas)
+template <class T>
+SpecialSquareMatrix<T>::SpecialSquareMatrix(SpecialSquareMatrix<T> const & a,  bool share)  
+  : mType(a.mType) , mNrows(a.mNrows) , mElems(a.mElems, share)  , mInfo(NULL)
+{
+  if (a.mInfo)  mInfo = new DVList(*(a.mInfo)); 
+}
+
+//! Sets or change the matrix size, specifying the new number of rows
+template <class T>
+void SpecialSquareMatrix<T>::SetSize(sa_size_t rowSize)
+{
+  throw  ForbiddenError("SpecialSquareMatrix<T>::SetSize() Only sub-classes should be instanciated");
+}
+
+//! Return the object  as a standard (TMatrix<T>) object
+template <class T>
+TMatrix<T> SpecialSquareMatrix<T>::ConvertToStdMatrix() const
+{
+  throw  ForbiddenError("SpecialSquareMatrix<T>::ConvertToStdMatrix() Only sub-classes should be instanciated");
+}
+
+
+
+
+//! Clone (duplicate) the SpecialSquareMatrix<T>  matrix \b a
+template <class T>
+void SpecialSquareMatrix<T>::Clone(SpecialSquareMatrix<T> const & a)
+{
+  if ((mType != C_UndefinedMatrix) && (a.mType != mType)) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::Clone()- Different type matrices");
+  mType = a.mType;
+  mNrows = a.mNrows;
+  mElems.Clone(a.mElems);
+  if (a.mInfo)  mInfo = new DVList(*(a.mInfo)); 
+}
+
+//! Share the data with \b a 
+template <class T>
+void SpecialSquareMatrix<T>::Share(SpecialSquareMatrix<T> const & a)
+{
+  if ((mType != C_UndefinedMatrix) && (a.mType != mType)) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::Share()- Different type matrices");
+  mType = a.mType;
+  mNrows= a.mNrows; 
+  mElems.Share(a.mElems);
+  if (a.mInfo)  mInfo = new DVList(*(a.mInfo)); 
+}
+
+//! Clone or share the data. Data is shared if \b is flagged as temporary
+template <class T>
+void SpecialSquareMatrix<T>::CloneOrShare(SpecialSquareMatrix<T> const & a)
+{
+  if ((mType != C_UndefinedMatrix) && (a.mType != mType)) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::CloneOrShare()- Different type matrices");
+  mType = a.mType;
+  mNrows = a.mNrows;
+  mElems.CloneOrShare(a.mElems);
+  if (a.mInfo)  mInfo = new DVList(*(a.mInfo)); 
+}
+
+
+
+//! Set all matrix elements to the value \b a
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::SetCst(T a)
+{
+  if (mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::SetCst(T )- (this) not allocated !");
+  mElems = a;
+  return (*this);
+}
+
+
+//! Element by element copy method 
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::Set(SpecialSquareMatrix<T> const & a)
+{
+  if (a.mType != mType) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::Set(a)- Different type matrices");
+  if (a.mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::Set(a)- Object a not allocated !");  
+  if (mElems.Size() < 1) {
+    mNrows = a.mNrows;
+    mElems.CloneOrShare(a.mElems);
+    return (*this);
+  }
+  if (mNrows != a.mNrows)  
+    throw SzMismatchError("SpecialSquareMatrix<T>::Set(a)- NRows() != a.NRows() ! ");
+  mElems = a.mElems;
+  return (*this);
+}
+
+//! Set the matrix elements using the values from the sequence \b seq
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::SetSeq(Sequence const & seq)
+{
+  if (mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::SetSeq()- (this) not allocated !");
+  for(size_t k=0; k<mElems.Size(); k++)  seq(k).Convert(mElems(k));
+  return (*this);
+}
+
+//! Add a constant value (\b b ) to all matrix elements 
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::AddCst(T b)
+{
+  if (mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::AddCst(T )- (this) not allocated !");
+  mElems += b;
+  return (*this);
+}
+
+//! Subtract a constant value (\b b ) from all matrix elements. fginv==true -> elem = b-elem
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::SubCst(T b, bool fginv)
+{
+  if (mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::Subst(T )- (this) not allocated !");
+  if (fginv) 
+    mElems = mElems.Sub(b, true);
+  else   mElems -= b;
+  return (*this);
+}
+
+//! Multiply all matrix elements by a constant value \b b . 
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::MulCst(T b)
+{
+  if (mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::MulCst(T )- (this) not allocated !");
+  mElems *= b; 
+  return (*this);
+}
+
+//! Divide all matrix elements by a constant value \b b . if fginv==true elem = b / elem 
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::DivCst(T b, bool fginv)
+{
+  if (mElems.Size() < 1)
+    throw SzMismatchError("SpecialSquareMatrix<T>::DivCst(T )- (this) not allocated !");
+  if (fginv)  mElems = mElems.Div(b, true);
+  else mElems /= b;
+  return (*this);
+}
+
+//! Element by element addition : elem(l,m) = elem(l,m) + b.elem(l,m)
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::AddElt(SpecialSquareMatrix<T> const& b)
+{
+  if (b.mType != mType) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::AddElt() Different type matrices");
+  if ((mNrows < 1)||(mNrows != b.mNrows)||(mElems.Size() < 1)||(mElems.Size() != b.mElems.Size()))
+    throw SzMismatchError("SpecialSquareMatrix<T>::AddElt()- Not allocated or different sizes !");
+  mElems += b.mElems;
+  return (*this);
+}
+
+//! Element by element subtraction : elem(l,m) = elem(l,m) - b.elem(l,m)
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::SubElt(SpecialSquareMatrix<T> const& b)
+{
+  if (b.mType != mType) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::SubElt() Different type matrices");
+  if ((mNrows < 1)||(mNrows != b.mNrows)||(mElems.Size() < 1)||(mElems.Size() != b.mElems.Size()))
+    throw SzMismatchError("SpecialSquareMatrix<T>::SubElt()- Not allocated or different sizes !");
+  mElems -= b.mElems;
+  return (*this);
+}
+//! Element by element multiplication : elem(l,m) = elem(l,m) * b.elem(l,m)
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::MulElt(SpecialSquareMatrix<T> const& b)
+{
+  if (b.mType != mType) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::MulElt() Different type matrices");
+  if ((mNrows < 1)||(mNrows != b.mNrows)||(mElems.Size() < 1)||(mElems.Size() != b.mElems.Size()))
+    throw SzMismatchError("SpecialSquareMatrix<T>::MulElt()- Not allocated or different sizes !");
+  mElems *= b.mElems;
+  return (*this);
+}
+//! Element by element divison : elem(l,m) = elem(l,m) / b.elem(l,m)
+template <class T>
+SpecialSquareMatrix<T>& SpecialSquareMatrix<T>::DivElt(SpecialSquareMatrix<T> const& b)
+{
+  if (b.mType != mType) 
+    throw TypeMismatchExc("SpecialSquareMatrix<T>::DivElt() Different type matrices");
+  if ((mNrows < 1)||(mNrows != b.mNrows)||(mElems.Size() < 1)||(mElems.Size() != b.mElems.Size()))
+    throw SzMismatchError("SpecialSquareMatrix<T>::DivElt()- Not allocated or different sizes !");
+  mElems /= b.mElems;
+  return (*this);
+}
+
+
+//! Return the minimum and the maximum values of the SpecialSquare Matrix elements 
+/*!
+  This method generates an exception (\c MathExc) if called for complex matrices
+*/
+template <class T>
+void SpecialSquareMatrix<T>::MinMax(T& min, T& max) const 
+{
+  const T * pe = mElems.Begin();
+  min = pe[0];
+  max = pe[0];
+  for(size_t k=0; k<mElems.Size(); k++) {
+    if (pe[k]<min)  min = pe[k];
+    else if (pe[k]>max)  max = pe[k];      
+  }
+  return;
+}
+
+DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
+void SpecialSquareMatrix< complex<r_4> >::MinMax(complex<r_4>& min, complex<r_4>& max) const
+{
+  throw MathExc("SpecialSquareMatrix< complex<r_4> >::MinMax(...) - No order in complex");
+}
+DECL_TEMP_SPEC  /* equivalent a template <> , pour SGI-CC en particulier */
+void SpecialSquareMatrix< complex<r_8> >::MinMax(complex<r_8>& min, complex<r_8>& max) const
+{
+  throw MathExc("SpecialSquareMatrix< complex<r_4> >::MinMax(...) - No order in complex");
+}
+
+//! Short ASCII representatio of the object pushed to stream \b os 
+template <class T>
+ostream& SpecialSquareMatrix<T>::Show(ostream& os) const
+{
+  os << typeid(*this).name() << " < " << typeid(T).name() 
+     << " > NRow=" << mNrows << " NbElem<>0 : " << Size() << endl;
+  return os;
+}
+
+
+}   // namespace SOPHYA
+
+///////////////////////////////////////////////////////////////
+#ifdef __CXX_PRAGMA_TEMPLATES__
+#pragma define_template SpecialSquareMatrix<uint_1>
+#pragma define_template SpecialSquareMatrix<uint_2>
+#pragma define_template SpecialSquareMatrix<uint_4>
+#pragma define_template SpecialSquareMatrix<uint_8>
+#pragma define_template SpecialSquareMatrix<int_1>
+#pragma define_template SpecialSquareMatrix<int_2>
+#pragma define_template SpecialSquareMatrix<int_4>
+#pragma define_template SpecialSquareMatrix<int_8>
+#pragma define_template SpecialSquareMatrix<r_4>
+#pragma define_template SpecialSquareMatrix<r_8>
+#pragma define_template SpecialSquareMatrix< complex<r_4> > 
+#pragma define_template SpecialSquareMatrix< complex<r_8> > 
+#ifdef SO_LDBLE128
+#pragma define_template SpecialSquareMatrix<r_16>
+#pragma define_template SpecialSquareMatrix< complex<r_16> > 
+#endif
+#endif
+
+#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
+namespace SOPHYA {
+template class SpecialSquareMatrix<uint_1>;
+template class SpecialSquareMatrix<uint_2>;
+template class SpecialSquareMatrix<uint_4>;
+template class SpecialSquareMatrix<uint_8>;
+template class SpecialSquareMatrix<int_1>;
+template class SpecialSquareMatrix<int_2>;
+template class SpecialSquareMatrix<int_4>;
+template class SpecialSquareMatrix<int_8>;
+template class SpecialSquareMatrix<r_4>;
+template class SpecialSquareMatrix<r_8>;
+template class SpecialSquareMatrix< complex<r_4> >;
+template class SpecialSquareMatrix< complex<r_8> >;
+#ifdef SO_LDBLE128
+template class SpecialSquareMatrix<r_16>;
+template class SpecialSquareMatrix< complex<r_16> >;
+#endif
+}
+#endif
+
Index: trunk/SophyaLib/TArray/spesqmtx.h
===================================================================
--- trunk/SophyaLib/TArray/spesqmtx.h	(revision 3809)
+++ trunk/SophyaLib/TArray/spesqmtx.h	(revision 3809)
@@ -0,0 +1,268 @@
+// This may look like C code, but it is really -*- C++ -*-
+// This code is part of the SOPHYA library
+//  (C) Univ. Paris-Sud   (C) LAL-IN2P3/CNRS   (C) IRFU-CEA
+//  (C) R. Ansari, C.Magneville    2009  
+
+#ifndef SPESQMTX_H_SEEN
+#define SPESQMTX_H_SEEN
+
+#include "ndatablock.h"
+#include "pexceptions.h"
+#include "utilarr.h"
+#include "tmatrix.h"
+
+
+namespace SOPHYA {
+
+//--- enum definissant le type de matrice carre 
+//! enum definition for the differents special square matrix types (or flavours)
+enum SpSqMtxType {
+  C_UndefinedMatrix = 0,        //!< Base type 
+  C_DiagonalMatrix = 1,         //!< Diagonal Matrix 
+  C_SymmetricMatrix = 4,        //!< Symmetric Matrix
+  C_LowerTriangularMatrix = 8,  //!< Lower Triangular Matrix
+  C_UpperTriangularMatrix = 9   //!< Upper Triangular Matrix
+};
+
+//-- forward declaration for the PPF handler 
+template <class T> class FIO_SpecialSquareMatrix;
+  
+template <class T>
+class SpecialSquareMatrix : public AnyDataObj {
+public :
+
+
+// Constructor - size 0, SetSize() should be called before the object is used
+explicit SpecialSquareMatrix(SpSqMtxType typ=C_UndefinedMatrix);
+
+// Constructeur avec specification du nombre de lignes 
+explicit SpecialSquareMatrix(sa_size_t rowSize, SpSqMtxType typ);
+
+// Constructor de copy (avec partage de reference de NDataBlock) 
+SpecialSquareMatrix(SpecialSquareMatrix<T> const & a,  bool share=true);  
+
+// Allocation / changement de taille - Doit etre redefinie par les classes derivees
+virtual void SetSize(sa_size_t rowSize);
+
+//! Alias for SetSize() 
+inline void ReSize(sa_size_t rowSize) 
+{
+  return SetSize(rowSize);
+}
+
+//! Alias for SetSize() 
+inline void ReSizeRow(sa_size_t rowSize) 
+{
+  return SetSize(rowSize);
+}
+
+// Retourne l'objet sous forme d'une TMatrix<T> 
+virtual TMatrix<T> ConvertToStdMatrix() const;
+
+// Clone : duplique l'objet et ses donnees 
+virtual void Clone(SpecialSquareMatrix<T> const & a);
+
+// Share : recopie l'objet a avec partage des donnees des elements  
+virtual void Share(SpecialSquareMatrix<T> const & a);
+
+// Operation de Share() si objet temporaire, Clone() sinon
+virtual void CloneOrShare(SpecialSquareMatrix<T> const & a);
+
+//! Flags the associated NDataBlock as temporary if \b fgtmp==true
+inline void SetTemp(bool fgtmp=false) const  { mElems.SetTemp(fgtmp); } 
+
+//! Return the matrix type (Diagonal, Triangular ...) 
+inline SpSqMtxType MtxType() const { return mType; }
+//! Return number of rows
+inline  sa_size_t NRows() const {return (int_4)mNrows;}
+//! Return number of colums ( = NRows())
+inline  sa_size_t NCols() const {return (int_4)mNrows;}
+
+//! Return the total number of non zero elements 
+inline sa_size_t Size() const { return mElems.Size();}
+
+//! Return reference to datablock containg matrix elements NDataBlock
+inline NDataBlock<T>& DataBlock()   { return mElems; }
+//! Return reference to datablock containg matrix elements NDataBlock
+inline const NDataBlock<T>& DataBlock() const  { return mElems; }
+
+//! Return a reference to the associated DVList object 
+inline DVList& Info()
+  { if (mInfo==NULL)  mInfo = new DVList;  return (*mInfo); }
+
+// Pour mettre tous les elements a la valeur a 
+virtual SpecialSquareMatrix<T>& SetCst(T a);
+
+//! Alias for SetT(a) - Set all matrix elements to the value \b a
+inline SpecialSquareMatrix<T>& SetT(T a)
+  { return SetCst(a); }
+
+//  Egalite/Recopie element par element 
+virtual SpecialSquareMatrix<T>& Set(SpecialSquareMatrix<T> const & a);
+
+//  Remplissage avec une sequence de valeurs
+virtual SpecialSquareMatrix<T>& SetSeq(Sequence const & seq);
+
+// Ajout de la constante "b: a  tous les elements
+virtual SpecialSquareMatrix<T>& AddCst(T b);
+
+// Soustrait "b" a  tous les elements - fginv==true -> elem = b-elem
+virtual SpecialSquareMatrix<T>& SubCst(T b, bool fginv=false);
+
+// Multiplie "b" a tous les elements 
+virtual SpecialSquareMatrix<T>& MulCst(T b);
+
+// Divise tous les elements par "b". fginv==true elem = b / elem 
+virtual SpecialSquareMatrix<T>& DivCst(T b, bool fginv=false);
+
+// Addition element par element:  elem(l,c) = elem(l,c) + b.elem(l,c)
+virtual SpecialSquareMatrix<T>& AddElt(SpecialSquareMatrix<T> const& b);
+
+// Soustraction element par element : elem(l,c) = elem(l,c) - b.elem(l,c)
+virtual SpecialSquareMatrix<T>& SubElt(SpecialSquareMatrix<T> const& b);
+
+// Multiplication element par element : elem(l,c) = elem(l,c) * b.elem(l,c)
+virtual SpecialSquareMatrix<T>& MulElt(SpecialSquareMatrix<T> const& b);
+
+// Division element par element : elem(l,c) = elem(l,c) / b.elem(l,c)
+virtual SpecialSquareMatrix<T>& DivElt(SpecialSquareMatrix<T> const& b);
+
+//--- Operateurs = (T b) , = (SpecialSquareMatrix<T> b)
+//! operator = a , to set all elements to the value \b a
+inline SpecialSquareMatrix<T>& operator = (T a)
+  {  return SetCst(a);  }
+//! operator = SpecialSquareMatrix<T> a , element by element copy operator
+inline SpecialSquareMatrix<T>& operator = (SpecialSquareMatrix<T> const & a)
+  {  return Set(a); }
+//! operator = Sequence seq , 
+inline SpecialSquareMatrix<T>& operator = (Sequence const & seq)
+  {  return SetSeq(seq); }
+
+//--- Operateur += -= *= /= (T b)
+//! operator += b , Add a constant value \b b to all elements
+inline SpecialSquareMatrix<T>& operator += (T b)
+  { return AddCst(b); }
+//! operator -= b , subtract a constant value \b b from all elements
+inline SpecialSquareMatrix<T>& operator -= (T b)
+  { return SubCst(b); }
+//! operator *= b , multiply all elements by  a constant value \b b 
+inline SpecialSquareMatrix<T>& operator *= (T b)
+  { return MulCst(b); }
+//! operator /= b , divide all elements by  a constant value \b b 
+inline SpecialSquareMatrix<T>& operator /= (T b)
+  { return DivCst(b); }
+
+//! Element access operator (R/W): access to the element at position \b index - No bound checking
+inline T& operator[](sa_size_t index) 
+  { return  mElems(index); }
+//! Element access operator (RO): access to the element at position \b index - No bound checking
+inline T operator[](sa_size_t index) const
+  { return  mElems(index); }
+
+// Valeur min et max des elements (sauf tableaux complexes -> exception)
+virtual void MinMax(T& min, T& max) const ;
+
+// Impression / representation texte sur le flot os
+virtual ostream& Show(ostream& os) const;
+
+friend class FIO_SpecialSquareMatrix<T>;
+//..............
+protected: 
+SpSqMtxType mType;        //!< Type of the square matrix (diagonal, symmetric, triangular)
+sa_size_t mNrows;         //!< size of the square matrix
+NDataBlock<T> mElems;     //!< Data block
+DVList* mInfo;            //!< Associated DVList structure
+};
+
+//! print (ascii representation) operator
+template <class T>
+inline ostream& operator << (ostream& os, const SpecialSquareMatrix<T>& a)
+                            { a.Show(os);    return(os);    }
+
+////////////////////////////////////////////////////////////////
+// Surcharge d'operateurs A (+,-,*,/) (T) x
+
+/*! \ingroup TArray \fn operator+(const SpecialSquareMatrix<T>&,T)
+  \brief Operator SpecialSquareMatrix = SpecialSquareMatrix + constant */
+template <class T> inline SpecialSquareMatrix<T> operator + (const SpecialSquareMatrix<T>& a, T b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.AddCst(b); return result; }
+
+/*! \ingroup TArray \fn operator+(T,const SpecialSquareMatrix<T>&)
+  \brief Operator SpecialSquareMatrix = constant + SpecialSquareMatrix */
+template <class T> inline SpecialSquareMatrix<T> operator + (T b,const SpecialSquareMatrix<T>& a)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.AddCst(b); return result; }
+
+/*! \ingroup TArray \fn operator-(const SpecialSquareMatrix<T>&,T)
+  \brief Operator SpecialSquareMatrix = SpecialSquareMatrix - constant */
+template <class T> inline SpecialSquareMatrix<T> operator - (const SpecialSquareMatrix<T>& a, T b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.SubCst(b); return result; }
+
+/*! \ingroup TArray \fn operator-(T,const SpecialSquareMatrix<T>&)
+  \brief Operator SpecialSquareMatrix = constant - SpecialSquareMatrix */
+template <class T> inline SpecialSquareMatrix<T> operator - (T b,const SpecialSquareMatrix<T>& a)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.SubCst(b, true); return result; }
+
+/*! \ingroup TArray \fn operator*(const SpecialSquareMatrix<T>&,T)
+  \brief Operator SpecialSquareMatrix = SpecialSquareMatrix * constant */
+template <class T> inline SpecialSquareMatrix<T> operator * (const SpecialSquareMatrix<T>& a, T b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.MulCst(b); return result; }
+
+/*! \ingroup TArray \fn operator*(T,const SpecialSquareMatrix<T>&)
+  \brief Operator SpecialSquareMatrix = constant * SpecialSquareMatrix */
+template <class T> inline SpecialSquareMatrix<T> operator * (T b,const SpecialSquareMatrix<T>& a)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.MulCst(b); return result; }
+
+/*! \ingroup TArray \fn operator/(const SpecialSquareMatrix<T>&,T)
+  \brief Operator SpecialSquareMatrix = SpecialSquareMatrix / constant */
+template <class T> inline SpecialSquareMatrix<T> operator / (const SpecialSquareMatrix<T>& a, T b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.DivCst(b); return result; }
+
+/*! \ingroup TArray \fn operator/(T,const SpecialSquareMatrix<T>&)
+  \brief Operator SpecialSquareMatrix = constant / SpecialSquareMatrix  */
+template <class T> inline SpecialSquareMatrix<T> operator / (T b, const SpecialSquareMatrix<T>& a)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.DivCst(b, true); return result; }
+
+
+////////////////////////////////////////////////////////////////
+// Surcharge d'operateurs C = A (+,-,&&,/) B
+
+/*! \ingroup TArray \fn operator+(const SpecialSquareMatrix<T>&,const SpecialSquareMatrix<T>&)
+  \brief Operator SpecialSquareMatrix = SpecialSquareMatrix + SpecialSquareMatrix (element by element addition) */
+template <class T>
+inline SpecialSquareMatrix<T> operator + (const SpecialSquareMatrix<T>& a,const SpecialSquareMatrix<T>& b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.AddElt(b); return result; }
+
+/*! \ingroup TArray \fn operator-(const SpecialSquareMatrix<T>&,const SpecialSquareMatrix<T>&)
+  \brief Operator SpecialSquareMatrix = SpecialSquareMatrix - SpecialSquareMatrix (element by element subtraction) */
+template <class T>
+inline SpecialSquareMatrix<T> operator - (const SpecialSquareMatrix<T>& a,const SpecialSquareMatrix<T>& b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.SubElt(b); return result; }
+
+/*! \ingroup TArray \fn operator && (const SpecialSquareMatrix<T>&,const SpecialSquareMatrix<T>&)
+  \brief Element by element multiplication of two arrays SpecialSquareMatrix = SpecialSquareMatrix * SpecialSquareMatrix */
+
+template <class T>
+inline SpecialSquareMatrix<T> operator && (const SpecialSquareMatrix<T>& a,const SpecialSquareMatrix<T>& b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.MulElt(b); return result; }
+
+/*! \ingroup TArray \fn operator / (const SpecialSquareMatrix<T>&,const SpecialSquareMatrix<T>&)
+  \brief Element by element division of two arrays SpecialSquareMatrix = SpecialSquareMatrix / SpecialSquareMatrix */
+template <class T>
+inline SpecialSquareMatrix<T> operator / (const SpecialSquareMatrix<T>& a,const SpecialSquareMatrix<T>& b)
+  { SpecialSquareMatrix<T> result; result.CloneOrShare(a); 
+    result.SetTemp(true);  result.DivElt(b); return result; }
+  
+}   // namespace SOPHYA
+
+#endif
Index: trunk/SophyaLib/TArray/spesqmtx_tsnl.h
===================================================================
--- trunk/SophyaLib/TArray/spesqmtx_tsnl.h	(revision 3809)
+++ trunk/SophyaLib/TArray/spesqmtx_tsnl.h	(revision 3809)
@@ -0,0 +1,28 @@
+#ifdef TSNMLUPG4
+/* This is a specific declaration file to cope with two level 
+   name look-up scheme , introduced by gcc >= 3.4 
+   THIS FILE SHOULD NOT NORMALLY BE INCLUDED IN USER PROGRAMS 
+*/
+
+using SOPHYA::SpecialSquareMatrix<T>::mType  ;
+using SOPHYA::SpecialSquareMatrix<T>::mNrows  ;
+using SOPHYA::SpecialSquareMatrix<T>::mElems  ;
+using SOPHYA::SpecialSquareMatrix<T>::mInfo  ;
+
+using SOPHYA::SpecialSquareMatrix<T>::SetSize  ;
+using SOPHYA::SpecialSquareMatrix<T>::ReSize  ;
+using SOPHYA::SpecialSquareMatrix<T>::ReSizeRow  ;
+
+using SOPHYA::SpecialSquareMatrix<T>::Info  ;
+
+using SOPHYA::SpecialSquareMatrix<T>::MtxType  ;
+using SOPHYA::SpecialSquareMatrix<T>::NRows  ;
+using SOPHYA::SpecialSquareMatrix<T>::NCols  ;
+using SOPHYA::SpecialSquareMatrix<T>::Size  ;
+using SOPHYA::SpecialSquareMatrix<T>::SetSeq  ;
+
+/*
+using SOPHYA::SpecialSquareMatrix<T> ;
+*/
+
+#endif  /* TSNMLUPG4 */
Index: trunk/SophyaLib/TArray/symmtx.h
===================================================================
--- trunk/SophyaLib/TArray/symmtx.h	(revision 3809)
+++ trunk/SophyaLib/TArray/symmtx.h	(revision 3809)
@@ -0,0 +1,226 @@
+// This may look like C code, but it is really -*- C++ -*-
+
+#ifndef SYMMTX_H_SEEN
+#define SYMMTX_H_SEEN
+
+#include "spesqmtx.h"
+
+// doit etre mis en dehors du namespace
+/*!
+  \class SOPHYA::SymmetricMatrix
+  \ingroup TArray
+  \brief Class representing a symmetric matrix. 
+
+  The symmetric matrix is represented in memory as column packed,
+  corresponding to the lower triangular part, as illustrated below for a 5x5 matrix.
+  \verbatim
+  5x5 symmetric.Matrix, Size= 5*(5+1)/2 = 15 independent elements (0 ... 14) 
+  | 0                  |
+  | 1   5              |
+  | 2   6   9          |
+  | 3   7   10  12     |
+  | 4   8   11  13  14 |
+  \endverbatim 
+
+  This class offers similar functionalities to the TArray<T> / TMatrix<T> classes, like 
+  reference sharing and counting, arithmetic operators ... However, this class has no 
+  sub matrix extraction method.
+*/
+
+namespace SOPHYA {
+  
+//! Class for inferior triangular matrix (base class for the class Alm)
+template <class T>
+class SymmetricMatrix : public SpecialSquareMatrix<T> {
+public :
+
+#include "spesqmtx_tsnl.h"
+
+//! Default constructor - TriangMatrix of size 0, SetSize() should be called before the object is used
+explicit SymmetricMatrix()
+  : SpecialSquareMatrix<T>(C_SymmetricMatrix)
+{
+
+}
+
+//! Instanciate a triangular matrix from the number of rows (rowSize must be > 0)
+explicit SymmetricMatrix(sa_size_t rowSize)
+  : SpecialSquareMatrix<T>(rowSize, C_SymmetricMatrix)
+{
+  if (rowSize < 1) 
+    throw ParmError("SymmetricMatrix<T>::SymmetricMatrix(rsz) rsz <= 0");
+  mElems.ReSize((rowSize*(rowSize+1)/2) ); 
+  mInfo = NULL;
+}
+
+//! Copy constructor (possibility of sharing datas)
+SymmetricMatrix(SymmetricMatrix<T> const & a,  bool share=false)  
+  : SpecialSquareMatrix<T>(a, share) 
+{
+}
+
+//! Copy constructor (possibility of sharing datas)
+SymmetricMatrix(SpecialSquareMatrix<T> const & a,  bool share=false)  
+  : SpecialSquareMatrix<T>(a, share) 
+{
+  if (a.MtxType() != C_SymmetricMatrix) 
+    throw TypeMismatchExc("SymmetricMatrix(a) a NOT a SymmetricMatrix");
+}
+
+/*!
+ \brief Create a lower triangular matrix from a square matrix. 
+ Elements above the diagonal are ignored.
+*/
+explicit SymmetricMatrix(TMatrix<T> const & mx)  
+  : SpecialSquareMatrix<T>(C_SymmetricMatrix)
+{
+  if ((mx.NRows() != mx.NCols()) || (mx.NRows() < 1)) 
+    throw ParmError("SymmetricMatrix<T>::(TMatrix<T> const & mx) mx not allocated OR NOT a square matrix");
+  SetSize(mx.NRows());
+  for(sa_size_t l=0; l<NRows(); l++) 
+    for(sa_size_t m=0; m<=l; m++) (*this)(l,m) = mx(l,m);
+}
+
+//! Sets or change the triangular matrix size, specifying the new number of rows
+virtual void SetSize(sa_size_t rowSize) 
+{
+  if (rowSize < 1) 
+    throw ParmError("SymmetricMatrix<T>::SetSize(rsz) rsz <= 0");
+  if (rowSize == mNrows)  return;
+  mNrows=rowSize;
+  mElems.ReSize(mNrows*(mNrows+1)/2);  
+}
+
+//! Return number of rows (for compatibility with the old TriangularMatrix interface)
+inline  sa_size_t rowNumber() const {return (int_4)mNrows;}
+
+//! Return the object (triangular matrix) as a standard square matrix
+virtual TMatrix<T> ConvertToStdMatrix() const
+{
+  if (mNrows < 1)
+    throw SzMismatchError("SymmetricMatrix<T>::ConvertToStdMatrix() (this) not allocated !");
+  SOPHYA::TMatrix<T> mx(NRows(), NRows());
+  for(sa_size_t l=0; l<NRows(); l++) 
+    for(sa_size_t m=0; m<=l; m++) mx(l,m) = mx(m,l) = (*this)(l,m);
+  return mx;
+}
+ 
+//--- Operateurs = (T b) , = (SymmetricMatrix<T> b)
+//! operator = a , to set all elements to the value \b a
+inline SymmetricMatrix<T>& operator = (T a)
+  {  SetCst(a);  return (*this);  }
+//! operator = SymmetricMatrix<T> a , element by element copy operator
+inline SymmetricMatrix<T>& operator = (SymmetricMatrix<T> const & a)
+  {  Set(a); return (*this); }
+//! operator = Sequence seq  
+inline SymmetricMatrix<T>& operator = (Sequence const & seq)
+  {  SetSeq(seq); return (*this); }
+
+
+//--- Operateurs d'acces aux elements 
+//! Element access operator (R/W): access to elements row \b r and column \b c
+inline T& operator()(sa_size_t r, sa_size_t c) 
+{ 
+  if ((r<0)||(r>=mNrows))  
+    throw RangeCheckError("DiagonalMatrix<T>::operator()(r,c) (r<0)||(r>=NRows())");
+  if (c>r) { sa_size_t rc = r; r=c; c=rc; }
+  // the inferior triangular part of the matrix is stored column by column
+  return(mElems(r+ mNrows*c-c*(c+1)/2));
+}
+//! Element access operator (RO): access to elements row \b l and column \b m
+inline T operator()(sa_size_t r, sa_size_t c) const
+{ 
+  if ((r<0)||(r>=mNrows))  
+    throw RangeCheckError("DiagonalMatrix<T>::operator()(r,c) (r<0)||(r>=NRows())");
+  if (c>r) { sa_size_t rc = r; r=c; c=rc; }
+  // the inferior triangular part of the matrix is stored column by column
+  return(mElems(r+ mNrows*c-c*(c+1)/2));
+}
+
+//! Return the pointer to the first non zero element in column \b j = &(tmmtx(j,j)) 
+inline const T* columnData(sa_size_t j)  const {return mElems.Begin()+(mNrows*j-j*(j-1)/2) ;}
+
+//! Return the pointer to the first non zero element in column \b j = &(tmmtx(j,j)) 
+inline T* columnData(sa_size_t j) {return mElems.Begin()+(mNrows*j-j*(j-1)/2) ;}
+
+//! compute the position of the element \b tm(i,j) relative to the first element
+inline sa_size_t indexOfElement(sa_size_t i,sa_size_t j) const 
+{
+  //  return(i*(i+1)/2+j);
+  // the (inferior triangular )matrix is stored column by column
+  return(i+ mNrows*j-j*(j+1)/2);
+}
+
+//! Triangular Matrix product (multiplication) : ret_matrix = (*this) * tmx
+TMatrix<T> Multiply(SymmetricMatrix<T> const & tmx) const
+{
+  if (NRows() != tmx.NRows()) 
+    throw SzMismatchError("Matrix<T>::Multiply(SymmetricMatrix<T> tmx): different sizes");
+// codage peu efficace : on utilise la multiplication de matrices generales ...
+  TMatrix<T> a = ConvertToStdMatrix();
+  TMatrix<T> b = tmx.ConvertToStdMatrix();
+  return (a.Multiply(b));
+}
+
+//! Matrix product (multiplication) : ret_matrix = (*this) * mx
+TMatrix<T> MultiplySG(TMatrix<T> const & mx) const 
+{
+  if (NCols() != mx.NRows())
+    throw SzMismatchError("SymmetricMatrix<T>::MultiplySG(TMatrix<T> mx): NCols()!=mx.NRows()");
+  TMatrix<T> a = ConvertToStdMatrix();
+  return a.Multiply(mx); 
+}
+
+//! Matrix product (multiplication) : ret_matrix = mx * (*this)
+TMatrix<T> MultiplyGS(TMatrix<T> const & mx) const 
+{
+  if (NRows() != mx.NCols())
+    throw SzMismatchError("SymmetricMatrix<T>::MultiplyGS(TMatrix<T> mx): NRows()!=mx.NCols()");
+  TMatrix<T> a = ConvertToStdMatrix();
+  return mx.Multiply(a); 
+}
+
+//! ASCII dump/print of the triangular matrix object (set nbLignes=-1 for dumping the complete matrix)
+ostream& Print(ostream& os, sa_size_t nbLignes=0) const
+{
+  os << "SymmetricMatrix< " << typeid(T).name() 
+     << " > NRow=" << mNrows << " NbElem<>0 : " << Size() << endl;
+  if (nbLignes == 0) return os;
+  if (nbLignes < 0 ) nbLignes = mNrows;
+  if (nbLignes > mNrows ) nbLignes = mNrows;
+  for (sa_size_t r=0; r<nbLignes; r++)  {
+    os << "Row[" << r << "]: " ;
+    for (sa_size_t c=0; c<NRows(); c++) 
+      os << " " << (*this)(r,c);
+    os << endl;
+  }
+  if (nbLignes < mNrows)  os << " ... ... ... " << endl;
+  return os;
+}
+
+protected: 
+};
+
+//----- Surcharge d'operateurs C = A * B (multiplication matricielle)
+/*! \ingroup TArray \fn operator*(const SymmetricMatrix<T>&,const SymmetricMatrix<T>&)
+  \brief * : SymmetricMatrix multiplication \b a and \b b */
+template <class T> 
+inline TMatrix<T> operator * (const SymmetricMatrix<T>& a, const SymmetricMatrix<T>& b)
+   { return(a.Multiply(b)); }
+
+/*! \ingroup TArray \fn operator*(const SymmetricMatrix<T>&,const TMatrix<T>&)
+  \brief * : Matrix multiplication SymmetricMatrix (\b a ) *  TMatrix<T> ( \b b ) */
+template <class T> 
+inline TMatrix<T> operator * (const SymmetricMatrix<T>& a, const TMatrix<T>& b)
+   { return(a.MultiplySG(b)); }
+
+/*! \ingroup TArray \fn operator*(const TMatrix<T>&,const SymmetricMatrix<T>&)
+  \brief * : Matrix multiplication TMatrix (\b a ) *  SymmetricMatrix<T> ( \b b ) */
+template <class T> 
+inline TMatrix<T> operator * (const TMatrix<T>& a, const SymmetricMatrix<T>& b)
+   { return(b.MultiplyGS(a)); }
+
+  
+}   // namespace SOPHYA
+
+#endif
Index: trunk/SophyaLib/TArray/tarrinit.cc
===================================================================
--- trunk/SophyaLib/TArray/tarrinit.cc	(revision 3791)
+++ trunk/SophyaLib/TArray/tarrinit.cc	(revision 3809)
@@ -5,4 +5,5 @@
 
 #include "fioarr.h"
+#include "fiospsqmtx.h"
 
 /*!
@@ -23,5 +24,9 @@
 //   V= 2.1 , Avr/Jul 2007: Constr.copie/extrac sous-tableaux thread-safe
 //   V= 2.2 , Fev/Nars 2010: Tableaux en long double 
-#define MOD_VERS   2.2
+//   V= 2.3 , Juillet 2010: Ajout des matrices carres Diagonal, Triangulaire, Symmetrique
+#define MOD_VERS   2.3
+
+// Fonction d'enregistrement des handlers matrices carres speciales 
+static void RegisterSSQM_PPFHandlers();
 
 TArrayInitiator::TArrayInitiator()
@@ -108,4 +113,6 @@
 #endif
 
+  RegisterSSQM_PPFHandlers();   // Enregistrement des gestionnaires PPF des matrices carrees speciales
+
   SophyaInitiator::RegisterModule("TArray", MOD_VERS);  // Module name and version number registration
 }
@@ -116,4 +123,85 @@
 }
 
+// --- Fonction d'initialisation / enregistrement des handlers PPF pour matrices carres speciales 
+static void RegisterSSQM_PPFHandlers()
+{
+// Enregistrement des gestionnaires PPF 
+  PPRegister(FIO_SpecialSquareMatrix<uint_1>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_1>, DiagonalMatrix<uint_1>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_1>, LowerTriangularMatrix<uint_1>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_1>, SymmetricMatrix<uint_1>);
+
+  PPRegister(FIO_SpecialSquareMatrix<uint_2>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_2>, DiagonalMatrix<uint_2>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_2>, LowerTriangularMatrix<uint_2>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_2>, SymmetricMatrix<uint_2>);
+
+  PPRegister(FIO_SpecialSquareMatrix<uint_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_4>, DiagonalMatrix<uint_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_4>, LowerTriangularMatrix<uint_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_4>, SymmetricMatrix<uint_4>);
+
+  PPRegister(FIO_SpecialSquareMatrix<uint_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_8>, DiagonalMatrix<uint_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_8>, LowerTriangularMatrix<uint_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<uint_8>, SymmetricMatrix<uint_8>);
+
+  PPRegister(FIO_SpecialSquareMatrix<int_1>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_1>, DiagonalMatrix<int_1>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_1>, LowerTriangularMatrix<int_1>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_1>, SymmetricMatrix<int_1>);
+
+  PPRegister(FIO_SpecialSquareMatrix<int_2>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_2>, DiagonalMatrix<int_2>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_2>, LowerTriangularMatrix<int_2>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_2>, SymmetricMatrix<int_2>);
+
+  PPRegister(FIO_SpecialSquareMatrix<int_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_4>, DiagonalMatrix<int_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_4>, LowerTriangularMatrix<int_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_4>, SymmetricMatrix<int_4>);
+
+  PPRegister(FIO_SpecialSquareMatrix<int_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_8>, DiagonalMatrix<int_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_8>, LowerTriangularMatrix<int_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<int_8>, SymmetricMatrix<int_8>);
+
+  PPRegister(FIO_SpecialSquareMatrix<r_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_4>, DiagonalMatrix<r_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_4>, LowerTriangularMatrix<r_4>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_4>, SymmetricMatrix<r_4>);
+
+  PPRegister(FIO_SpecialSquareMatrix<r_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_8>, DiagonalMatrix<r_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_8>, LowerTriangularMatrix<r_8>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_8>, SymmetricMatrix<r_8>);
+
+  PPRegister(FIO_SpecialSquareMatrix< complex<r_4> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_4> >, DiagonalMatrix< complex<r_4> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_4> >, LowerTriangularMatrix< complex<r_4> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_4> >, SymmetricMatrix< complex<r_4> >);
+
+  PPRegister(FIO_SpecialSquareMatrix< complex<r_8> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_8> >, DiagonalMatrix< complex<r_8> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_8> >, LowerTriangularMatrix< complex<r_8> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_4> >, SymmetricMatrix< complex<r_8> >);
+
+#ifdef SO_LDBLE128
+  PPRegister(FIO_SpecialSquareMatrix<r_16>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_16>, DiagonalMatrix<r_16>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_16>, LowerTriangularMatrix<r_16>);
+  DObjRegister(FIO_SpecialSquareMatrix<r_16>, SymmetricMatrix<r_16>);
+
+  PPRegister(FIO_SpecialSquareMatrix< complex<r_16> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_16> >, DiagonalMatrix< complex<r_16> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_16> >, LowerTriangularMatrix< complex<r_16> >);
+  DObjRegister(FIO_SpecialSquareMatrix< complex<r_16> >, SymmetricMatrix< complex<r_16> >);
+#endif 
+
+  return;
+}
+
+
+
 
 // On met un objet initiator en statique, pour les loaders qui savent 
Index: trunk/SophyaLib/TArray/triangmtx.h
===================================================================
--- trunk/SophyaLib/TArray/triangmtx.h	(revision 3791)
+++ 	(revision )
@@ -1,202 +1,0 @@
-// This may look like C code, but it is really -*- C++ -*-
-
-#ifndef TRIANGMTX_H_SEEN
-#define TRIANGMTX_H_SEEN
-
-#include <typeinfo>
-#include "ndatablock.h"
-#include "pexceptions.h"
-
-// doit etre mis en dehors du namespace
-/*!
-  \class SOPHYA::TriangularMatrix
-  \ingroup TArray
-  \brief Class for inferior triangular matrix (base class for the class Alm)
-  The inferior triangular matrix is represented in memory as column packed,
-  as illustrated below for a 5x5 triangular matrix.
-  \verbatim
-  5x5 Inf.Triang.Matrix, Size= 15 elements (0 ... 14) 
-  | 0                  |
-  | 1   5              |
-  | 2   6   9          |
-  | 3   7   10  12     |
-  | 4   8   11  13  14 |
-  \endverbatim 
-*/
-
-namespace SOPHYA {
-  
-//! Class for inferior triangular matrix (base class for the class Alm)
-template <class T>
-class TriangularMatrix {
-public :
-
-//! Default constructor
-TriangularMatrix()
-  : long_diag_(0)
-{
-}
-
-//! instanciate a triangular matrix from the number of rows
-TriangularMatrix(sa_size_t rowSize)
-  : long_diag_(rowSize) 
-{
-  elem_.ReSize((rowSize*(rowSize+1)/2) ); 
-}
-
-//! Copy constructor (possibility of sharing datas)
-TriangularMatrix(const TriangularMatrix<T>& a,  bool share=false)  
-  : long_diag_(a.long_diag_) , elem_(a.elem_, share)  
-{
-}
-
-//! resize the matrix with a new number of rows
-inline void ReSizeRow(sa_size_t rowSize) 
-{
-  long_diag_=(uint_4)rowSize;
-  elem_.ReSize(long_diag_*(long_diag_+1)/2);
-}
-
-TriangularMatrix<T>& SetT(T a)
-  {
-    if (long_diag_ < 1)
-    throw RangeCheckError("TriangularMatrix<T>::SetT(T )  - TriangularMatrix not dimensionned ! ");
-    elem_ = a;
-    return (*this);
-  }
-
-//! () operator : access to elements row \b l and column \b m
-inline T& operator()(sa_size_t l, sa_size_t m) 
-  {
-      return  elem_(indexOfElement(l,m));
-  }
-
-inline T& operator()(sa_size_t index) 
-  {
-      return  elem_(index);
-  }
-
-
-//! () operator : access to elements row \b l and column \b m
-inline T const& operator()(sa_size_t l, sa_size_t m) const 
-  {
-      return *(elem_.Begin()+ indexOfElement(l,m));
-  }
-
-inline T const& operator()(sa_size_t index) const 
-  {
-      return *(elem_.Begin()+ index);
-  }
-
-TriangularMatrix<T>& Set(const TriangularMatrix<T>& a)
-  {
-    if (this != &a) 
-      {
-	if (a.Size() < 1) 
-	  throw RangeCheckError(" TriangularMatrix<T>::Set()- Array a not allocated ! ");
-      }
-    if (Size() < 1)  CloneOrShare(a);
-    else CopyElt(a);
-    return(*this);
-  }
-
-inline TriangularMatrix<T>& operator = (const TriangularMatrix<T>& a)
-                                                           {return Set(a);}
-
-TriangularMatrix<T>& CopyElt(const  TriangularMatrix<T>& a)
-  {
-    if (Size() < 1) 
-      throw RangeCheckError("TriangularMatrix<T>::CopyElt(const TriangularMatrix<T>& )  - Not Allocated Array ! ");
-    if (Size() != a.Size() )
-      throw(SzMismatchError("TriangularMatrix<T>::CopyElt(const TriangularMatrix<T>&) SizeMismatch")) ;
-    long_diag_ = a.long_diag_;
-    sa_size_t k;
-    for (k=0; k< Size(); k++) elem_(k) = a.elem_(k);
-    return(*this);
-  }
-
-void CloneOrShare(const  TriangularMatrix<T>& a)
-  {
-    long_diag_ = a.long_diag_;
-    elem_.CloneOrShare(a.elem_);
-  }
-
-
-//! Return number of rows
-inline  sa_size_t rowNumber() const {return (int_4)long_diag_;}
-
-//! Return size of the total array
-  inline sa_size_t Size() const {return elem_.Size();}
-
-  inline bool CheckRelativeIndices(sa_size_t l, sa_size_t m) const 
-   {
-     if ( l < m ) 
-       {
-	 throw RangeCheckError("TriangularMatrix<T>::CheckRelativeIndices: indices out of range " );
-       }
-     return true;
-   }
-  inline bool CheckAbsoluteIndice(sa_size_t l, sa_size_t m) const 
-   {
-     if ( indexOfElement(l,m) >= elem_.Size() )
-       {
-	 throw RangeCheckError("TriangularMatrix<T>::CheckAbsoluteIndice: indices out of range " );
-       }
-   }
-  inline bool CheckAbsoluteIndice(sa_size_t ind) const 
-   {
-     if ( ind >= elem_.Size() )
-       {
-	 throw RangeCheckError("TriangularMatrix<T>::CheckAbsoluteIndice: indices out of range " );
-       }
-   }
-
-  //! ASCII dump of the matrix (set nbLignes=-1) for dumping the complete matrix
-void Print(ostream& os, sa_size_t nbLignes=0) const
-  {
-    os << "TriangularMatrix< " << typeid(T).name() 
-       << " > NRow=" << long_diag_ << " NbElem<>0 : " << Size() << endl;
-    if (nbLignes == 0) return;
-    if (nbLignes < 0 ) nbLignes = long_diag_;
-    if (nbLignes > long_diag_ ) nbLignes = long_diag_;
-    for (sa_size_t k=0; k < nbLignes; k++)  {
-      os << "L[" << k << "]: " ;
-      for (sa_size_t kc = 0; kc <= k ; kc++) 
-	os << " " << elem_(indexOfElement(k,kc));
-      os << endl;
-    }
-    if (nbLignes < long_diag_)  os << " ... ... ... " << endl;
-    return;
-  }
-
-inline void Print(sa_size_t nbLignes=0) const  { Print(cout, nbLignes); }
-
-
-//! Return the pointer to the first non zero element in column \b j = &(tmmtx(j,j)) 
-inline const T* columnData(sa_size_t j)  const {return elem_.Begin()+(long_diag_*j-j*(j-1)/2) ;}
-
-//! Return the pointer to the first non zero element in column \b j = &(tmmtx(j,j)) 
-inline T* columnData(sa_size_t j) {return elem_.Begin()+(long_diag_*j-j*(j-1)/2) ;}
-
-//! compute the address of an element in the single array representing the matrix
-inline sa_size_t indexOfElement(sa_size_t i,sa_size_t j) const 
-{
-  //  return(i*(i+1)/2+j);
-  // the (inferior triangular )matrix is stored column by column
-  return(i+ long_diag_*j-j*(j+1)/2);
-}
-
-private: 
-
-sa_size_t long_diag_;    //!< size of the square matrix
-NDataBlock<T> elem_;  //!< Data block
-
-  };
-
-template <class T>
-inline ostream& operator << (ostream& os, const TriangularMatrix<T>& a)
-                            { a.Print(os, 0);    return(os);    }
-  
-}   // namespace SOPHYA
-
-#endif
Index: trunk/SophyaLib/TArray/trngmtx.h
===================================================================
--- trunk/SophyaLib/TArray/trngmtx.h	(revision 3809)
+++ trunk/SophyaLib/TArray/trngmtx.h	(revision 3809)
@@ -0,0 +1,234 @@
+// This may look like C code, but it is really -*- C++ -*-
+
+#ifndef TRNGMTX_H_SEEN
+#define TRNGMTX_H_SEEN
+
+#include "spesqmtx.h"
+
+// doit etre mis en dehors du namespace
+/*!
+  \class SOPHYA::LowerTriangularMatrix
+  \ingroup TArray
+  \brief Class representing a lower (inferior) triangular matrix. This is the 
+  base class for the Alm<T> class (in module Samba), representing complex 
+  coefficients of spherical harmonic decomposition.
+
+  The lower triangular matrix is represented in memory as column packed,
+  as illustrated below for a 5x5 triangular matrix.
+  \verbatim
+  5x5 Inf.Triang.Matrix, Size= 5*(5+1)/2 = 15 elements (0 ... 14) 
+  | 0                  |
+  | 1   5              |
+  | 2   6   9          |
+  | 3   7   10  12     |
+  | 4   8   11  13  14 |
+  \endverbatim 
+
+  This class offers similar functionalities to the TArray<T> / TMatrix<T> classes, like 
+  reference sharing and counting, arithmetic operators ... However, this class has no 
+  sub matrix extraction method.
+*/
+
+namespace SOPHYA {
+  
+//! Class for inferior triangular matrix (base class for the class Alm)
+template <class T>
+class LowerTriangularMatrix : public SpecialSquareMatrix<T> {
+public :
+
+#include "spesqmtx_tsnl.h"
+
+//! Default constructor - TriangMatrix of size 0, SetSize() should be called before the object is used
+explicit LowerTriangularMatrix()
+  : SpecialSquareMatrix<T>(C_LowerTriangularMatrix)
+{
+  mZeros = T(0);
+}
+
+//! Instanciate a triangular matrix from the number of rows (rowSize must be > 0)
+explicit LowerTriangularMatrix(sa_size_t rowSize)
+  : SpecialSquareMatrix<T>(rowSize, C_LowerTriangularMatrix)
+{
+  if (rowSize < 1) 
+    throw ParmError("LowerTriangularMatrix<T>::LowerTriangularMatrix(rsz) rsz <= 0");
+  mElems.ReSize((rowSize*(rowSize+1)/2) ); 
+  mInfo = NULL;
+  mZeros = T(0);
+}
+
+//! Copy constructor (possibility of sharing datas)
+LowerTriangularMatrix(LowerTriangularMatrix<T> const & a,  bool share=false)  
+  : SpecialSquareMatrix<T>(a, share) 
+{
+}
+
+//! Copy constructor (possibility of sharing datas)
+LowerTriangularMatrix(SpecialSquareMatrix<T> const & a,  bool share=false)  
+  : SpecialSquareMatrix<T>(a, share) 
+{
+  if (a.MtxType() != C_LowerTriangularMatrix) 
+    throw TypeMismatchExc("LowerTriangularMatrix(a) a NOT a LowerTriangularMatrix");
+  mZeros = T(0);
+}
+
+/*!
+ \brief Create a lower triangular matrix from a square matrix. 
+ Elements above the diagonal are ignored.
+*/
+explicit LowerTriangularMatrix(TMatrix<T> const & mx)  
+  : SpecialSquareMatrix<T>(C_LowerTriangularMatrix)
+{
+  if ((mx.NRows() != mx.NCols()) || (mx.NRows() < 1)) 
+    throw ParmError("LowerTriangularMatrix<T>::(TMatrix<T> const & mx) mx not allocated OR NOT a square matrix");
+  SetSize(mx.NRows());
+  for(sa_size_t l=0; l<NRows(); l++) 
+    for(sa_size_t m=0; m<=l; m++) (*this)(l,m) = mx(l,m);
+  mZeros = T(0);
+}
+
+//! Sets or change the triangular matrix size, specifying the new number of rows
+virtual void SetSize(sa_size_t rowSize) 
+{
+  if (rowSize < 1) 
+    throw ParmError("LowerTriangularMatrix<T>::SetSize(rsz) rsz <= 0");
+  if (rowSize == mNrows)  return;
+  mNrows=rowSize;
+  mElems.ReSize(mNrows*(mNrows+1)/2);  
+}
+
+//! Return number of rows (for compatibility with the old TriangularMatrix interface)
+inline  sa_size_t rowNumber() const {return (int_4)mNrows;}
+
+//! Return the object (triangular matrix) as a standard square matrix
+virtual TMatrix<T> ConvertToStdMatrix() const
+{
+  if (mNrows < 1)
+    throw SzMismatchError("LowerTriangularMatrix<T>::ConvertToStdMatrix() (this) not allocated !");
+  SOPHYA::TMatrix<T> mx(NRows(), NRows());
+  for(sa_size_t l=0; l<NRows(); l++) 
+    for(sa_size_t m=0; m<=l; m++) mx(l,m) = (*this)(l,m);
+  return mx;
+}
+ 
+//--- Operateurs = (T b) , = (LowerTriangularMatrix<T> b)
+//! operator = a , to set all elements to the value \b a
+inline LowerTriangularMatrix<T>& operator = (T a)
+  {  SetCst(a);  return (*this);  }
+//! operator = LowerTriangularMatrix<T> a , element by element copy operator
+inline LowerTriangularMatrix<T>& operator = (LowerTriangularMatrix<T> const & a)
+  {  Set(a); return (*this); }
+//! operator = Sequence seq  
+inline LowerTriangularMatrix<T>& operator = (Sequence const & seq)
+  {  SetSeq(seq); return (*this); }
+
+
+//--- Operateurs d'acces aux elements 
+//! Element access operator (R/W): access to elements row \b r and column \b c
+inline T& operator()(sa_size_t r, sa_size_t c) 
+{ 
+  if ((r<0)||(r>=mNrows))  
+    throw RangeCheckError("DiagonalMatrix<T>::operator()(r,c) (r<0)||(r>=NRows())");
+  if (c>r) { mZeros = T(0); return mZeros; }
+  // the (inferior triangular )matrix is stored column by column
+  return(mElems(r+ mNrows*c-c*(c+1)/2));
+}
+//! Element access operator (RO): access to elements row \b l and column \b m
+inline T operator()(sa_size_t r, sa_size_t c) const
+{ 
+  if ((r<0)||(r>=mNrows))  
+    throw RangeCheckError("DiagonalMatrix<T>::operator()(r,c) (r<0)||(r>=NRows())");
+  if (c>r) { mZeros = T(0); return mZeros; }
+  // the (inferior triangular )matrix is stored column by column
+  return(mElems(r+ mNrows*c-c*(c+1)/2));
+}
+
+//! Return the pointer to the first non zero element in column \b j = &(tmmtx(j,j)) 
+inline const T* columnData(sa_size_t j)  const {return mElems.Begin()+(mNrows*j-j*(j-1)/2) ;}
+
+//! Return the pointer to the first non zero element in column \b j = &(tmmtx(j,j)) 
+inline T* columnData(sa_size_t j) {return mElems.Begin()+(mNrows*j-j*(j-1)/2) ;}
+
+//! compute the position of the element \b tm(i,j) relative to the first element
+inline sa_size_t indexOfElement(sa_size_t i,sa_size_t j) const 
+{
+  //  return(i*(i+1)/2+j);
+  // the (inferior triangular )matrix is stored column by column
+  return(i+ mNrows*j-j*(j+1)/2);
+}
+
+//! Triangular Matrix product (multiplication) : ret_matrix = (*this) * tmx
+LowerTriangularMatrix<T> Multiply(LowerTriangularMatrix<T> const & tmx) const
+{
+  if (NRows() != tmx.NRows()) 
+    throw SzMismatchError("Matrix<T>::Multiply(LowerTriangularMatrix<T> tmx): different sizes");
+// codage peu efficace : on utilise la multiplication de matrices generales ...
+  TMatrix<T> a = ConvertToStdMatrix();
+  TMatrix<T> b = tmx.ConvertToStdMatrix();
+  LowerTriangularMatrix<T> ret(a.Multiply(b));
+  ret.SetTemp(true);
+  return ret;
+}
+
+//! Matrix product (multiplication) : ret_matrix = (*this) * mx
+TMatrix<T> MultiplyTG(TMatrix<T> const & mx) const 
+{
+  if (NCols() != mx.NRows())
+    throw SzMismatchError("LowerTriangularMatrix<T>::MultiplyTG(TMatrix<T> mx): NCols()!=mx.NRows()");
+  TMatrix<T> a = ConvertToStdMatrix();
+  return a.Multiply(mx); 
+}
+
+//! Matrix product (multiplication) : ret_matrix = mx * (*this)
+TMatrix<T> MultiplyGT(TMatrix<T> const & mx) const 
+{
+  if (NRows() != mx.NCols())
+    throw SzMismatchError("LowerTriangularMatrix<T>::MultiplyGT(TMatrix<T> mx): NRows()!=mx.NCols()");
+  TMatrix<T> a = ConvertToStdMatrix();
+  return mx.Multiply(a); 
+}
+
+//! ASCII dump/print of the triangular matrix object (set nbLignes=-1 for dumping the complete matrix)
+ostream& Print(ostream& os, sa_size_t nbLignes=0) const
+{
+  os << "LowerTriangularMatrix< " << typeid(T).name() 
+     << " > NRow=" << mNrows << " NbElem<>0 : " << Size() << endl;
+  if (nbLignes == 0) return os;
+  if (nbLignes < 0 ) nbLignes = mNrows;
+  if (nbLignes > mNrows ) nbLignes = mNrows;
+  for (sa_size_t k=0; k < nbLignes; k++)  {
+    os << "L[" << k << "]: " ;
+    for (sa_size_t kc = 0; kc <= k ; kc++) 
+      os << " " << mElems(indexOfElement(k,kc));
+    os << endl;
+  }
+  if (nbLignes < mNrows)  os << " ... ... ... " << endl;
+  return os;
+}
+
+protected: 
+mutable T mZeros; 
+};
+
+//----- Surcharge d'operateurs C = A * B (multiplication matricielle)
+/*! \ingroup TArray \fn operator*(const LowerTriangularMatrix<T>&,const LowerTriangularMatrix<T>&)
+  \brief * : LowerTriangularMatrix multiplication \b a and \b b */
+template <class T> 
+inline LowerTriangularMatrix<T> operator * (const LowerTriangularMatrix<T>& a, const LowerTriangularMatrix<T>& b)
+   { return(a.Multiply(b)); }
+
+/*! \ingroup TArray \fn operator*(const LowerTriangularMatrix<T>&,const TMatrix<T>&)
+  \brief * : Matrix multiplication LowerTriangularMatrix (\b a ) *  TMatrix<T> ( \b b ) */
+template <class T> 
+inline TMatrix<T> operator * (const LowerTriangularMatrix<T>& a, const TMatrix<T>& b)
+   { return(a.MultiplyTG(b)); }
+
+/*! \ingroup TArray \fn operator*(const TMatrix<T>&,const LowerTriangularMatrix<T>&)
+  \brief * : Matrix multiplication TMatrix (\b a ) *  LowerTriangularMatrix<T> ( \b b ) */
+template <class T> 
+inline TMatrix<T> operator * (const TMatrix<T>& a, const LowerTriangularMatrix<T>& b)
+   { return(b.MultiplyGT(a)); }
+
+  
+}   // namespace SOPHYA
+
+#endif
