NSF grant

Linda Kaufman

Solving Symmetric Band Systems and Other Problems in Fiber Optic Design

Description for layman

The chemical composition of the layers of an optical fiber determines the optical properties of a fiber. Because a fiber suitable for underwater transmission, as an example, would not be appropriate for a local area network or to splice into an existing network to restore a degraded signal, there is not one fiber design that is suitable for all fibers. Our aim in this project is to produce tools for specifying the properties of a fiber and then to suggest several designs for the fiber using the latest models. Early on it was realized that especially for a model that considers fiber wrapped around a spool we needed to solve thousands of linear systems that were large, symmetric, and banded. One could use a classical technique that ignored symmetry, but using symmetry not only cuts down on the computational requirements but also could be used to refine the course of the sequence of linear systems, to hasten the convergence of the of the underlying nonlinear eigenvalue problem. Such a situation also exists in the design of the cavity of a linear collider and the design of bridges(Think of the TacomaNarrowsBridge collapse.), oil platforms, and buildings that might be affected by the environment. Part of this project involves investigating the properties and defining a good implementation for the algorithm for factoring a symmetric, banded system devised by the principal investigator.

The primary problem when using algorithms that preserve symmetry on a banded matrix is that they can require permutations of the matrix to limit element growth. These permutations mean that the banded matrix develops wings which can expand so that the band structure is essentially ruined. The PI has developed an algorithm for retracting the wings at each step and this project involves developed the algorithm into good software and studying its properties. Because this banded structure comes up often when approximating partial differential equations in various engineering problems, it is a rather useful tool.

Impact on other disciplines

Obviously this project has implications for designing optical fibers and structural engineering problems. When a signal is sent down a fiber, it tends to spread out over various frequencies- it tends to disperse. The dispersion can be computed using the second derivatives of the positive eigenvalues of an eigenvalue problem that may be nonlinear because of boundary conditions. To solve the inverse problem of finding the chemical composition means we need the gradient of the dispersion with respect to the design parameters, the composition and widths of the layers of the fiber, to feed into an optimizer. Because the dispersion is itself a derivative one can interchange the order of differentiating and we have found that interchanging decreases the number of systems that must be solved and allows for expanded use of the sparse structure of the Jacobian of the model matrix.

Students:

Serafim Stamatiswill be a junior in computer science and his job was to make the output of the testing code for the retraction algorithm more attractive. He had to learn the input/output system for C. From January 2007 to May 2007 he spent about 5 hours per week at $10/hour on the project.

Seok-Min was a South Korean sophomore transfer student, who probably has the highest mathematical aptitude among our computer science undergraduates. He was writing code for producing matrices to plug in to my algorithm for calculating the dispersion slope and its gradient with respect to design parameters for a particular fiber optics design. The dispersion slope uses the third derivative of the eigenvalues with respect to frequency defined by Maxwell's equations. His project involved finding the total derivatives of a multi variable function and modifying a code of 3000 lines that I had been working on that involved second derivatives. As he said, not only the length but the complexity of the code was a big step from the classroom programs. Seok-Min joined the project in January 2007 and will continue in the coming year.

During the school year he works about 10 hours a week exclusive of exam period and during the summer 20 hours/week at $10/ hour

Diamond is a Hispanic woman who will be a senior in Mathematics at WPU this fall. From January 2007 to May 2007 she spent about 5 hours per week at $10/hour on the project of computing the condition number of a symmetric banded matrix using the retraction algorithm developed by the PI and we hope this will become her senior thesis pending approval of the mathematics department.

Yomi is a black American and will be a senior this coming year at WPU. He joined the project in May 2007 for the summer and works about 5 hours per week at $10/hour. He has been writing code to compare my symmetric banded factorization with LAPACK's symmetric non banded and unsymmetric banded matrices on the symmetric indefinite matrices in U. of Florida sparse matrix collection.

Philip is an Asian American and is receiving a B.S. degree August 2007 in Computer Science at WPU. He joined the project in May 2007 and will leave in August. He has been taking symmetric indefinite banded matrices from the University of Florida sparse matrix collection, putting them through a reverse Cuthill-McGee algorithm to minimize their bandwidths and then comparing various versions of the retraction algorithm and other codes for factoring symmetric banded matrices. He has also been writing a special version of the retraction algorithm for 7 diagonal matrices for the fiber design problem. He has been working 20 hours/week at $10/ hour

Orielsy is an Hispanic veteran of the war in Iraq who was the top freshman in my beginning programming class in the fall of 2006 at WPU. He used F2C to translate some of the Fortran code I had written in 1999 at Bell Labs for finding eigenvalues and finding properties for a particular optical design. He is augmenting the code to specify the matrix dictated by Marcuse's model for a fiber wrapped around a spool and will link it into the code written by Phillip Kang and myself for 7 diagonal symmetric linear systems arising when solving the underlying nonlinear eigenvalue problem and determining the dispersion and its gradient for a particular design. Orielsy joined the project in January 2007 and will continue in the coming year.

During the school year he works about 10 hours a week exclusive of exam period and during the summer 20 hours/week at $10/ hour

Training and development:

For the 6 students who have worked on the project, the experience has been an eye opener. Only one of the students had ever worked on any program that was more than a toy class room project of say 100 lines. Here they were working on programs that had physical significance and were thousands of lines where each line depended on mathematics and translating that mathematics correctly. I tried to impress upon them that their codes had to be self documenting so they could be understood by other people and that correctness did not mean a lower grade but the collapse of an oil platform. They began to appreciatethe significance of the concepts in their math courses like singularity, total differentiation, and multivariable calculus. They learned how to get around some of the features of Visual C++ so they could handle large matrices and they learned how cache misses could affect the computation Some students learned more in a few months than from years of courses. They all developed a sense of pride that they could handle what was required.

I saw each student about an hour and a half a week. Some were independent and I was amazed; others needed a lot of attention and needed to have concepts explained again and again and I thought the time required with them was an order of magnitude more than if I done it myself.

Outreach activities

In April 2007 I gave a talk at WilliamPatersonUniversity on Research Day. Because the audience was primarily faculty members from other disciplines, I decided to jettison my slides on the details of the retraction algorithm.Itried to impress on them the utility of linear algebra in their fields and in their lives citing examples from latent semantic indexing, statistics, google, medical imaging, music, the design of optical fibers, bridge construction, etc. I mentioned that these applications led to linear systems like they had solved in their algebra classes but with thousands of variables and that one had to be worried then about the amount of computation time, which meant taking into consideration any structure the problem might possess. Moreover one had to be concerned about getting erroneous results because of the finite precision of the machine.

Project activities

There are two main objectives of RUI grant: (1) studying methods for solving the inverse problem to determine the chemical composition of a spooled optical fiber with particular optical properties which involves solving many banded symmetric linear systems and (2) analyzing and implementing various versions of the retraction algorithm devised in the Summer of 2005 by the Principal Investigator (PI) for factoring a symmetric band matrix which may be indefinite. Some of the students on the grant have been working on the fiber optic project making more efficient the packagewritten at Bell Labs and adding features to handle the dispersion slope and the fiber wrapped around the spool as in Marcuse’s model. Other students are comparingfor accuracy and timing various versions of the retraction code and snap- back from TelAvivUniversity on problems from sparse matrix collections. Results have been reported in seminars at Bell Labs, TelAvivUniversity, and at WilliamPatersonUniversity. A paper on calculating dispersion was published in a fiber optics journal and one on the algorithm for symmetric banded systems was published in a numerical linear algebra journal. A seminar at Fiber optics solution is anticipated.

Findings

(1)The construction of the matrices of the fiber optics problems can be sped up significantly by using the fact that the elements of the matrix can be written by a linear combination of separable terms where each term is a product of a function that depends only on frequency and another that depends only on the refractive index profile( the design parameters) of the optical fiber. This simplifies the calculation of derivatives with respect to frequency which are needed to calculation the dispersion- the spread of a signal along various wavelengths and the gradient of the various optical properties with respect to the design parameters.

The dispersion can be computed using the second derivatives of the positive eigenvalues of an eigenvalue problem that may be nonlinear because of boundary conditions. To solve the inverse problem of finding the chemical composition means we need the gradient of the dispersion with respect to the design parameters, the composition and widths of the layers of the fiber, to feed into an optimizer. Because the dispersion is itself a derivative one can interchange the order of differentiating and we have found that interchanging decreases the number of systems that must be solved and allows for expanded use of the sparse structure of the derivative of the main matrix derived from the model.

(2) One of the referees for the optics paper pointed out that designers were now using more complex models(even with complex arithmetic) and wondered if my formulae carried over to these full vector models. I realized that I should emphasize that my formulae for calculating the dispersion and its gradient were independent of the approximation technique and the complexity of the model. Moreover the new models meant that the bandwidths would be larger so that the retraction algorithmw ould have a wider market.

(3) In the retraction algorithm- if the matrix is symmetric but indefinite, whenever there is a 2 x 2 pivot( because a 1x 1 pivot would cause element growth) pivoting for stability may cause the matrix to develop wings, which the retraction algorithms eliminates with a series of transformations which could be orthogonal or could be stabilized elementary. Worst case analysis suggests that using orthogonal transformations could increase the computation time by almost a factor of 3 but conventional wisdom claims it would be more stable. In practice we have found that often there is little increase in the computation time and that the element growth is often greater with the orthogonal transformation and the relative residuals have not been uniformly smaller. It seems that with orthogonal transformations especially on random examples there has been a decrease in the number of 2 x 2’s and the wings have been smaller.

(4) Our experimentation indicates that the time for the retraction algorithm is significantly more than the time for blocked unsymmetric banded DGBTRF from LAPACK when the bandwidth becomes greater than say 500 on our machines. For our algorithm half the space is used for computation and half to store the information when the wings appear. Dividing this space into 2 separate matrices decreases cache problems and has lead to a 50% decrease in time. For unsymmetric band LU. 2/3 the space is used for computation and 1/3 the space is used for storage for the forward solve.

In the problem is mostly positive definite, one could create a blocking scheme, but it is the presence of the wings that prevents a full blocked algorithm, We have tried to use level 2 BLAS when possible- especially when using elementary transformations during the retraction phase but that produces only about a 10% change. However our experimentation has lead us to consider a frontal scheme attacking multiple wings at once.

(5) In our current version, unlike LAPACK’s symmetric indefinite non banded code, pivoting is not done during a 1 x 1 pivot. We have implemented a program that mimics the protocol in LAPACK for the banded case. It requires 50% more space and can be about almost twice as expensive, but the lower bound on element growth does not seem in practice to be important.

(6)

Human resources

I have been particular mindful of the student/minority aspect of this RUI grant. One of the students is a black American from Newark and two are Hispanic. I had particularly chosen younger students because I knew they would still be around and continuation was important.Moreover, I have found that when students receive encouragement, whether monetary or just relating to them as responsible adults, they tend to rise to the occasion.

The students who have worked for me the previous summers with support from the university have given talks in the department and poster sessions at St. Joseph’s student symposium. I will encourage the students under the grant to give talks this winter.