Mathematica程序 计算两线段之间的最短距离

Mathematica程序计算两线段之间的最短距离(*输入四个点坐标*)(*将坐标输成7/2形式可以得到准确解,输成3.5只能得到有一定精确度的数值解.*)(*如果两线段分别为AB,CD,按照以下规则输入坐标:{xA,yA,xB,yB,xC,yC,xD,yD}或者{{xA,yA},{xB,yB},{xC,yC},{xD,yD}}或者两种混和输入,必须保证按照顺序输入*)zuobiao = {{0, 0}, {13, 9}, {6, 6}, {0, 8}};zuobiao = Partition[Flatten[zuobiao], 2];(*下面进行坐标旋转,目的是让两条直线均不与坐标轴平行,省去以后很多分情况的麻烦,一劳永逸*)xuanzj = {0, \[Pi]/6, \[Pi]/4};brr = {};err = {};juli = {};frr = {};chuizuf =.;For[i = 1, i <= Length[xuanzj], i++,For[j = 1, j <= Length[zuobiao], j++,brr = Append[brr, {zuobiao[[j]][[1]]*Cos[xuanzj[[i]]] +zuobiao[[j]][[2]]*Sin[xuanzj[[i]]], -zuobiao[[j]][[1]]*Sin[xuanzj[[i]]] + zuobiao[[j]][[2]]*Cos[xuanzj[[i]]]}];];If[brr[[1]][[1]] != brr[[2]][[1]] &&brr[[3]][[1]] != brr[[4]][[1]] && brr[[1]][[2]] != brr[[2]][[2]] && brr[[3]][[2]] != brr[[4]][[2]],\[Alpha] = -xuanzj[[i]];Break[],brr = {}];];(*下面判断两直线是否平行*)pingxing =.;If[(zuobiao[[1]][[1]] - zuobiao[[2]][[1]])*(zuobiao[[3]][[2]] -zuobiao[[4]][[2]]) == (zuobiao[[1]][[2]] -zuobiao[[2]][[2]])*(zuobiao[[3]][[1]] - zuobiao[[4]][[1]]),If[(zuobiao[[1]][[1]] - zuobiao[[3]][[1]])*(zuobiao[[3]][[2]] -zuobiao[[4]][[2]]) == (zuobiao[[1]][[2]] -zuobiao[[3]][[2]])*(zuobiao[[3]][[1]] - zuobiao[[4]][[1]]), pingxing = 2(*直线重合*),pingxing = 1(*直线平行*)],pingxing = 0(*直线相交*)];(*下面是两直线重合时,先判断两条线段是否有重合部分,再计算距离*)If[pingxing == 2,drr = Transpose[brr];If[Partition[Sort[drr[[1]], Less],2] == {Sort[{drr[[1]][[1]], drr[[1]][[2]]}, Less],Sort[{drr[[1]][[3]], drr[[1]][[4]]}, Less]} ||Partition[Sort[drr[[1]], Less],2] == {Sort[{drr[[1]][[3]], drr[[1]][[4]]}, Less],Sort[{drr[[1]][[1]], drr[[1]][[2]]}, Less]},juli = Simplify[Norm[{Sort[drr[[1]], Less][[2]] - Sort[drr[[1]], Less][[3]],Sort[drr[[2]], Less][[2]] - Sort[drr[[2]], Less][[3]]}]];frr = {{{Sort[drr[[1]], Less][[2]],Sort[drr[[2]], Less][[2]]}, {Sort[drr[[1]], Less][[3]],Sort[drr[[2]], Less][[3]]}}};Print["准确解"];Print[juli];Print["数值解"];Print[juli // N],Print["两线段有一部分重合,距离为0"]]](*下面是两直线相交时,先计算两直线交点坐标,再判断两条线段是否相交*)(*下面定义一个函数jiaodian,用来计算两直线交点坐标*)jiaodian[zuobiao_] :=Module[{arr},arr = Flatten[zuobiao];{-(-arr[[3]] arr[[5]] arr[[2]] + arr[[3]] arr[[7]] arr[[2]] +arr[[1]] arr[[5]] arr[[4]] - arr[[1]] arr[[7]] arr[[4]] +arr[[1]] arr[[7]] arr[[6]] - arr[[3]] arr[[7]] arr[[6]] -arr[[1]] arr[[5]] arr[[8]] +arr[[3]] arr[[5]] arr[[8]])/(arr[[5]] arr[[2]] -arr[[7]] arr[[2]] - arr[[5]] arr[[4]] + arr[[7]] arr[[4]] -arr[[1]] arr[[6]] + arr[[3]] arr[[6]] + arr[[1]] arr[[8]] -arr[[3]] arr[[8]]), -(-arr[[3]] arr[[2]] arr[[6]] +arr[[7]] arr[[2]] arr[[6]] + arr[[1]] arr[[4]] arr[[6]] -arr[[7]] arr[[4]] arr[[6]] + arr[[3]] arr[[2]] arr[[8]] -arr[[5]] arr[[2]] arr[[8]] - arr[[1]] arr[[4]] arr[[8]] +arr[[5]] arr[[4]] arr[[8]])/(arr[[5]] arr[[2]] -arr[[7]] arr[[2]] - arr[[5]] arr[[4]] + arr[[7]] arr[[4]] -arr[[1]] arr[[6]] + arr[[3]] arr[[6]] + arr[[1]] arr[[8]] -arr[[3]] arr[[8]])}];(*下面是两直线相交时,先计算两直线交点坐标,再判断两线段是否相交,如果相交,输出"两线段相交,最短距离为0",否则继续接着的步骤*) xiangjiao =.;If[pingxing == 0,brr = Append[brr, jiaodian[brr]];(*计算两直线交点坐标*)If[(brr[[1]][[1]] - brr[[5]][[1]])*(brr[[5]][[1]] - brr[[2]][[1]]) >=0 && (brr[[3]][[1]] - brr[[5]][[1]])*(brr[[5]][[1]] -brr[[4]][[1]]) >= 0,Print["两线段相交,最短距离为0"];xiangjiao = 1,xiangjiao = 0];];(*下面先定义一个函数chuizu,用来计算点到直线的垂足坐标*) chuizu[zuobiao_] :=Module[{arr},arr = Flatten[zuobiao];{-(-arr[[1]] arr[[3]]^2 + 2 arr[[1]] arr[[3]] arr[[5]] -arr[[1]] arr[[5]]^2 - arr[[3]] arr[[2]] arr[[4]] +arr[[5]] arr[[2]] arr[[4]] - arr[[5]] arr[[4]]^2 +arr[[3]] arr[[2]] arr[[6]] - arr[[5]] arr[[2]] arr[[6]] +arr[[3]] arr[[4]] arr[[6]] + arr[[5]] arr[[4]] arr[[6]] -arr[[3]] arr[[6]]^2)/(arr[[3]]^2 - 2 arr[[3]] arr[[5]] +arr[[5]]^2 + arr[[4]]^2 - 2 arr[[4]] arr[[6]] +arr[[6]]^2), -(-arr[[1]] arr[[3]] arr[[4]] +arr[[1]] arr[[5]] arr[[4]] + arr[[3]] arr[[5]] arr[[4]] -arr[[5]]^2 arr[[4]] - arr[[2]] arr[[4]]^2 +arr[[1]] arr[[3]] arr[[6]] - arr[[3]]^2 arr[[6]] -arr[[1]] arr[[5]] arr[[6]] + arr[[3]] arr[[5]] arr[[6]] +2 arr[[2]] arr[[4]] arr[[6]] -arr[[2]] arr[[6]]^2)/(arr[[3]]^2 - 2 arr[[3]] arr[[5]] +arr[[5]]^2 + arr[[4]]^2 - 2 arr[[4]] arr[[6]] + arr[[6]]^2)} ];(*下面是两线段线不相交时,计算各个端点到所对线段的最短距离*) crr = {{brr[[1]], brr[[3]], brr[[4]]}, {brr[[2]], brr[[3]],brr[[4]]}, {brr[[3]], brr[[1]], brr[[2]]}, {brr[[4]], brr[[1]],brr[[2]]}};If[xiangjiao == 0,For[i = 1, i <= 4, i++,chuizuf = chuizu[crr[[i]]];If[(crr[[i]][[2]][[1]] - chuizuf[[1]])/(crr[[i]][[2]][[1]] -crr[[i]][[3]][[1]]) >=0 && (crr[[i]][[3]][[1]] - chuizuf[[1]])/(crr[[i]][[3]][[1]] -crr[[i]][[2]][[1]]) >= 0,frr = Append[frr, {crr[[i]][[1]], chuizuf}];juli =Append[juli, Simplify[Norm[crr[[i]][[1]] - chuizuf]]](*A的垂足F 在CD内,取线段AF距离*)];If[(crr[[i]][[2]][[1]] - chuizuf[[1]])/(crr[[i]][[2]][[1]] -crr[[i]][[3]][[1]]) <0 && (crr[[i]][[3]][[1]] - chuizuf[[1]])/(crr[[i]][[3]][[1]] -crr[[i]][[2]][[1]]) >= 0,frr = Append[frr, {crr[[i]][[1]], crr[[i]][[2]]}];juli =Append[juli,Simplify[Norm[crr[[i]][[1]] - crr[[i]][[2]]]]](*A的垂足F在C的外侧,取线段AC距离*)];If[(crr[[i]][[2]][[1]] - chuizuf[[1]])/(crr[[i]][[2]][[1]] -crr[[i]][[3]][[1]]) >=0 && (crr[[i]][[3]][[1]] - chuizuf[[1]])/(crr[[i]][[3]][[1]] -crr[[i]][[2]][[1]]) < 0,frr = Append[frr, {crr[[i]][[1]], crr[[i]][[3]]}];juli =Append[juli,Simplify[Norm[crr[[i]][[1]] - crr[[i]][[3]]]]](*A的垂足F在D的外侧,取线段AD距离*)];];Print["准确解"];Print[Sort[juli, Less][[1]]];Print["数值解"];Print[Sort[juli // N][[1]]];];(*下面是作图*)For[j = 1, j <= Length[brr], j++,err = Append[err, {brr[[j]][[1]]*Cos[\[Alpha]] +brr[[j]][[2]]*Sin[\[Alpha]], -brr[[j]][[1]]*Sin[\[Alpha]] +brr[[j]][[2]]*Cos[\[Alpha]]}];];err = Partition[err, 2];If[frr != {},grr = Table[Norm[frr[[i]][[1]] - frr[[i]][[2]]] // N, {i, 1, Length[frr]}];For[i = 1, i <= Length[grr], i++,If[grr[[i]] < grr[[1]],grr[[1]] = grr[[i]];frr[[1]] = frr[[i]]]];err = Append[err, frr[[1]]];];Print["下图中蓝色为已知的两条直线,红色为连接两条已知线段的最短线段"];ListLinePlot[err, AspectRatio -> Automatic,PlotStyle -> {Blue, Blue, Red}, AxesOrigin -> {0, 0},PlotRange -> All]。