Inverse Problems solved with a Calculus-level (Programming) Language

By Phil B Brubaker

Optimal Designs Enterprise

Licensed according to this deed.

Published on

Abstract

Inverse Problems

“You know what you want, just don’t know how to get there.”

Inverse Problems (IPs) are an important special set of problems that fit the statement “You know what you want, you just don’t know how to get there.”  A proper directive can get an Inverse Problem solved in hours!  Here are some examples of Inverse Problems:

  1. Law enforcement: Shells found at scene where did it come from?
  2. Airplane crash with wreckage all over the place.  How did these parts get where they lay?
  3. Missile target: Have target, how to get missile there?
  4. Want a ‘black box’ to have an efficiency of 54.321%.  How to design/build such a black box?
  5. Car seat storage H x W x D slot, how to design seat so as it will fit into slot while maximizing seat comfort?

A Parameter Estimation for an Inverse Problem is solved using the Calculus-level Find statement shown here:

Find a   ooo   To Match Error

If a Bound Value Problem then: Find a, ydot0, y2dot0   ooo   To Match Error

Where ‘a’ may be a vector with ‘n’ parts, a1, a2, a3,…an;
ydot0, y2dot0,etc. are derivatives at independent variable = 0; and,
‘error’ is the objective function.

If the IP problem contains any Ordinary or Partial Differential Equations (ODEs or PDEs) , the ‘find’ statement is wrapped around an integrate & integration statement in order to solve the ODEs or PDEs while finding the best ‘a’ parameter(s) for the given problem.

The ‘a’ parameter(s) are varied to fit one’s ‘m’ data points that make up the objective function, error.  This technique can vary as many parameters as you want; e.g. 5 or 50 or 50,000.  If there are less equations than parameters m < n, this would be classified as an under-determined system of equations.  If there are more equations than parameters, m > n, this would be an over-determined system.  Under- or Over-determined systems might force one to switch solvers to do the job.

Drug Development
from the pharmacokinetics arena

Problem Description

A drug company is trying to develop a drug with a given/desired half-life.  They know what they want, they just don’t know how to get there.  A true inverse problem!  Drug companies have several basic components that need to be combined together to determine what combination will provide the best fit to their desired half-life.  Here is the basic code for finding the best combination.

Computer App for Solving these Problems

A Calculus-level (Programming) Language is required to solve this type of problem.  Please download the free FortranCalculus compiler, install, and execute some Demo files in order to get an idea of the process that one must do to solve your problem.  If you want the manual, download the Basic Manual.

The following is just a rough sketch of necessary code to solve such a problem

Computer Code

The ‘Find’ statement is the work horse of a Calculus-level compiler.  It calls ones math model as many times as necessary in order to converge on a solution.  It varies your parameters, in this case (vol1, vol2, & maybe vol3), as it calls your math model.  The ‘in’ phase tells the name of math model routine.  The ‘by’ phase tells what solver to use, Asolver here; e.g. could be AJAX.  And the ‘to’ phase tells what the objective function is; ‘match’ means all following variables must equal zero, ‘g’ variable in this case.

FIND vol1, vol2,vol3;      IN drugs;       BY Asolver;       TO MATCH g

Global all
Problem halflife
  prnt = 0:   desired = 1.234        ! desired half-life time
  R = 0.9998            ! .02% Improvement factor
  nQty = 2:   nDrugs = ???
  vol1 = 1:   vol2 = 0:   vol3 = 0
  Do 20 i=1, nDrugs
    Do 20 j = i+1, nDrugs
      Find vol1,vol2; in drugs;  by Ajax(cntl);  to Match g
      If(g .lt. R*errmin) then
        Errmin=g: idrug=i:  jdrug=j
      endif
20 continue 
  nQty = 3:   vol1 = 1:   vol2 = 0:   vol3 = 0
  if( nDrugs .lt. nQty) goto 99
  Do 30 i=1, nDrugs
    Do 30 j = i+1, nDrugs
      Do 30 k = j+1, nDrugs
      Find vol1,vol2,vol3; in drugs;  by Asolver(cntl);  to Match g
      If(g .lt. R*errmin) then
        Errmin=g: idrug=i:  jdrug=j:  kdrug =k
      endif
30 continue 
99 continue 
nQty = 2:   prnt = 1   ! Print final results
  i = idrug:  j = jdrug:  k = kdrug:   vol1=1:  vol2=0:  vol3=0
  if( kdrug .lt. 1) then
      Find vol1,vol2; in drugs;  by Asolver(cntl);  to Match g
  else
    nQty = 3
      Find vol1,vol2,vol3; in drugs;  by Ajax(cntl);  to Match g
  endif
end
controller cntl for ‘Asolver’
   summary = prnt    ! solver only prints when summary=1
end
model drugs
  dimension cl( ’nDrugs’)  ! ‘nDrugs’ characteristic go here
  data cl/ 0.1234, 987.54, etc. /  ! ‘cl’ = clearance factor
 
  p1 = vol1 / cl(i):   p2 = vol2 / cl(j):   p3 = 0
  if( nQty .eq. 3)  p3 = vol3 / cl(k)
  halfLif = .693 * (p1 + p2 + p3)
  g = (desired – halfLif)**2
  cost = ???:     sellPric = ???:    profit = sellPric - cost
end
 

Findings

This drug program would execute all ‘nDrugs’ models and determine which model was best.  No output during the first pass.  Once the best was determined, the ‘prnt’ switch would be turned on to allow for a summary output to be printed showing the resulting parameters vol1, vol2, & maybe vol3.

Practical mode: round solution vales to practical values; e.g. say vol2 = 12.345678… but from a practical value you want vol1 = 12.4.

nQty = 3:   prnt = 1   ! Print final results
! i = idrug:  j = jdrug:  k = kdrug       ! set these index findings on next line.
i = 33:  j = 22:  k = 55:   vol2=12.4:  vol1=0:  vol3=0
 
      Find vol1,vol3; in drugs;  by Ajax(cntl);  to Match g

You may need to do this rounding again with vol1 and/or vol3.

Your Turn!

Have an Inverse Problem to solve?  Please state it here and use graphs, pictures, etc. to get your problem well stated and understood by those reading it.

Cite this work

Researchers should cite this work as follows:

  • Phil B Brubaker (2021), "Inverse Problems solved with a Calculus-level (Programming) Language," https://simiode.org/resources/8359.

    BibTex | EndNote

Tags