Plotting a curve from the intersection points of a family of two curves

I have two state equations.

f1 = 1000. (12 - (1 - u) x2);
f2 = 1000. ((1 - u) x1 - 8/100 x2 );


I did two different stream plots for state-variables x2 and x1 on the x and y axis respectively. The first stream plot was based on u = 0 and the second one was based on u = 1. Then I superimposed the stream plots using the command Show[{plt0, plt1}].

 plt0 = StreamPlot[Flatten[{f1, f2} /. u -> 0], {x2, 0, 30}, {x1, 0, 10},
StreamPoints -> 500, StreamScale -> Large,
StreamStyle -> {Blue, Thickness[0.001], Dotted}, Frame -> True,
FrameLabel -> {Row[{Subscript["x", 2], " - Capacitor  Voltage "}],
Row[{Subscript["x", 1], " - Inductor  Voltage "}]},
GridLines -> Automatic, GridLinesStyle -> Directive[Gray, Dashed],
BaseStyle -> {"Times", 24}, ImageSize -> 500]


 plt1 = StreamPlot[Flatten[{f1, f2} /. u -> 1], {x2, 0, 30}, {x1, 0, 10},
StreamPoints -> 500, StreamScale -> Large,
StreamStyle -> {Red, Thickness[0.001], Dotted}, Frame -> True,
FrameLabel -> {Row[{Subscript["x", 2], " - Capacitor  Voltage "}],
Row[{Subscript["x", 1], " - Inductor  Voltage "}]},
GridLines -> Automatic, GridLinesStyle -> Directive[Gray, Dashed],
BaseStyle -> {"Times", 24}, ImageSize -> 500]


    Show[{plt0, plt1}]


Is it possible to draw a curve at the intersection points of the two superimposed streamplots and how?

• While asking a question it is always better to supply some code snippet that can reproduce or pinpoint the problem you are facing. This spares others from spending time to reformulate a problem for which you may already have a code at hand. Dec 19, 2012 at 0:34
• I'm not really understand your question. What is the curve you want? Or you just want an arbitrary curve which goes through some intersection point in your given plot? Dec 19, 2012 at 3:08
• Virtually every point is the intersection of two stream lines of the vector field at u=0 and u=1. StreamPlot only gives a sample of all the stream lines. So a single curve won't go through all the intersection points (of all the stream lines). Are you looking for something like a wave front? As I understand StreamPlot, the actual stream lines drawn do not have a regular pattern. A line stops when an area gets too crowded, and vice versa, one starts when an area is too empty. The intersections of two StreamPlot's are not very telling. Dec 19, 2012 at 3:17
• Actually, for this particualar system of state equations, each line from the first plot uniquely intersects with one and only one line in the second plot. So, by superimposing them, i wanted to draw curve through those unique intersecting points. Any help would be appreaciated. I was superimposing them through this command Show[{plt0, plt1}]
– ram
Dec 19, 2012 at 3:39
• I'm sorry I don't understand. Looking at the lines in the combined plot, I see each blue line crossing several red lines, not a unique one. Dec 19, 2012 at 12:50

2 Answers

As @MichaelE2 noted in the comments, the lines produced by StreamPlot are just "a sample from all the stream lines" and the number of them and the extent of each depend on the option settings as well as finer implementation details. Yet, OP seems to be interested in the intersection of lines in a particular instance of StreamPlot output with specific option settings. The following is an attempt to get the intersection points of the lines shown in a graph produced by StreamPlot by

1. getting the coordinates of Line primitives,
2. producing an interpolating function for each line taking into account the range of its x coordinates,
3. solving for the first intersection point for each pair of lines from two different stream plots, and
4. adding the points thus obtained to the original stream plots using Epilog.

First, the final result (using StreamPoints -> 10):

Some preps:

  (* options to be used later *)
strmPltOpts = {Frame -> True,
FrameLabel -> {Row[{Subscript["x", 2], " - Capacitor Voltage "}],
Row[{Subscript["x", 1], " - Inductor Voltage "}]},
GridLines -> Automatic, GridLinesStyle -> Directive[Gray, Dashed],
BaseStyle -> {"Times", 24}, ImageSize -> 500};
pltOpts = {PlotRange -> {0, 10}, Frame -> True,
FrameLabel -> {Row[{Subscript["x", 2], " - Capacitor  Voltage "}],
Row[{Subscript["x", 1], " - Inductor  Voltage "}]},
GridLines -> Automatic, PlotStyle -> Thick,
GridLinesStyle -> Directive[Gray, Dashed],
BaseStyle -> {"Times", 24}, ImageSize -> 500, AspectRatio -> 1};


Few helper functions to get the lines, roots, and intersections:

  ClearAll[getRootsF, getStrmPltLinesF, getIntrsctnsF];
getStrmPltLinesF[plt_, x_] :=  Cases[Normal[plt[[1]]],
Line[z_] :>  ConditionalExpression[Interpolation[z][x],
Min[First@Transpose[z]] <= x <= Max[First@Transpose[z]]], Infinity];
getRootsF :=  Quiet[If[# ===  Undefined, {},
(soln = FindRoot[ First@#,
{x, Sequence @@ {First@#, First@#, Last@#} &@Last[#]}];
{{x /. soln, Last[First[#]] /. soln}, {x /. soln, First[First[#]] /. soln}})]] &;
getIntrsctnsF := First[Flatten[Quiet[getRootsF /@ #] /. {} -> Sequence[] /.
{{x_, y_}, {x_, z_}} /; Chop[y + z] != 0 :>
Sequence[] /. {_, _?Negative} :> Sequence[], {2}]] &;


Step 0: StreamPlot both original plots with the option setting StreamScale->None:

 plt0B = StreamPlot[Flatten[{f1, f2} /. u -> 0], {x2, 0, 30}, {x1, 0, 10},
StreamPoints -> 10, StreamScale -> None,
StreamStyle -> {Blue, Thickness[0.001], Dotted}, Evaluate@strmPltOpts];
plt1B = StreamPlot[Flatten[{f1, f2} /. u -> 1], {x2, 0, 30}, {x1, 0, 10},
StreamPoints -> 10, StreamScale -> None,
StreamStyle -> {Red, Thickness[0.001], Dotted},  Evaluate@strmPltOpts];
Show[{plt0B, plt1B}]


Steps 1 and 2: Use the helper function getStrmPltLinesF[strmplt, x] to get an interpolating function for each line in strmplt. Here is how these interpolating functions look when used in Plot:

 plt0Ba = Quiet@ Plot[Evaluate@getStrmPltLinesF[plt0B, x], {x, 0, 30},
ColorFunction -> "LakeColors", Evaluate[pltOpts]];
plt1Bb = Quiet@ Plot[Evaluate@getStrmPltLinesF[plt1B, x], {x, 0, 30},
ColorFunction -> "SolarColors", Evaluate[pltOpts]];
Show[{plt0Ba, plt1Bb}]


Step 3: get the points of intersection for every pair of lines using the helper function getIntrsctnsF:

 pairwiseDiffs = Join @@ Outer[Subtract, getStrmPltLinesF[plt0B, x],
getStrmPltLinesF[plt1B, x]];
intrsctnPoints = getIntrsctnsF@pairwiseDiffs;


Last step: change the options settings back to the original and add the intersection points using Epilog to one of the two original stream plots:

 plt0Bfinal = StreamPlot[Flatten[{f1, f2} /. u -> 0], {x2, 0, 30}, {x1, 0, 10},
StreamPoints -> 10, StreamScale -> Large,
StreamStyle -> {Blue, Thickness[0.001], Dotted},  Evaluate@strmPltOpts,
Epilog -> {Green, PointSize[.02], Point[intrsctnPoints]}];
plt1Bfinal = StreamPlot[Flatten[{f1, f2} /. u -> 1], {x2, 0, 30}, {x1, 0, 10},
StreamPoints -> 10, StreamScale -> Large,
StreamStyle -> {Red, Thickness[0.001], Dotted}, Evaluate@strmPltOpts];
Show[{plt0Bfinal, plt1Bfinal}]


to get the first picture above.

• +1 for the effort. By the way, welcome to the 20K club! :-) Dec 19, 2012 at 21:49
• @Mr.W thank you for the vote and the welcome.
– kglr
Dec 19, 2012 at 21:59
• Thanks a lot Mr Wizard. I'll take a look. I cam work with this right now.
– ram
Dec 19, 2012 at 22:19

This is intended as a complement to kguler's response. I add to those stream plots the curves where intersections are parallel and orthogonal respectively.

x[u_, x1_, x2_] := 1000* (12 - (1 - u) x2)
y[u_, x1_, x2_] := 1000*((1 - u) x1 - 8/100 x2)


Compute the implicit forms of the parallel and perpendicularly intersecting point sets (both are curves).

perppoly = {x[0, x1, x2], y[0, x1, x2]}.{x[1, x1, x2],
y[1, x1, x2]}

(* 12000000 (12 - x2) - 80000 (x1 - (2 x2)/25) x2 *)

parapoly =
First[GroebnerBasis[
Thread[{x[0, x1, x2], y[0, x1, x2]} -
lambda*{x[1, x1, x2], y[1, x1, x2]}], lambda]]

(* -150 x1 + x2^2 *)


Create the plots.

cperp = ContourPlot[{x[0, x1, x2], y[0, x1, x2]}.{x[1, x1, x2],
y[1, x1, x2]} == 0, {x2, 0, 30}, {x1, 0, 10},
ContourStyle -> {Thick, Green, Dashed}];

cparal = ContourPlot[parapoly == 0, {x2, 0, 30}, {x1, 0, 10},
ContourStyle -> {Thick, Black, Dotted}];


Now show the whole thing.

Show[plt0, plt1, cperp, cparal]


• Thanks a ton Mr Wizard and Daniel and all others who helped out. I think the black dotted curve is what i was looking for. this Surely helps. And i sincerely apologize if my questions were confusing but my understanding of mathematica is very basic and so i probably couldn't clarify enough.
– ram
Dec 20, 2012 at 8:11